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EXECUTIVE  SUMMARY 


Thrust  vectoring  is  one  means  of  improving  the  maneuverability  and  take-off  and 
landing  performance  of  aircraft.  However,  thrust  vectoring,  especially  in  close  proximity  to 
the  ground,  can  increase  significantly  the  complexity  of  the  flow  field  around  the  aircraft, 
and  make  aircraft  design  more  difficult.  To  achieve  good  aircraft  designs  in  spite  of  the 
complex  flow,  it  is  necessary  to  have  an  accurate  flow  prediction  method.  Methods  based 
on  the  Reynolds-averaged  Navier-Stokes  equations  can  make  these  predictions,  if  the 
numerical  discretization  and  the  turbulence  modeling  are  sufficiently  accurate.  This  work 
addresses  the  issue  of  turbulence  modeling  for  the  inclined  impinging  jets  associated  with 
the  use  of  thrust  reversing  in  ground  effect. 

Large  eddy  simulations  (LES)  and  physical  experiments  were  conducted  to  obtain 
mean  flow  and  turbulence  data.  These  data  were  analyzed  to  identify  the  dominant 
physical  processes  which  govern  turbulence  in  these  flows.  Streamline  curvature  and 
vortex  stretching  were  identified  as  the  two  most  significant  influences  on  turbulence  in 
these  flows  which  are  not  treated  adequately  in  the  standard  k-e  model.  Modifications 
which  account  for  curvature  and  vortex  stretching  were  added  to  the  standard  k-£  model. 
The  new  coefficients  in  the  resulting  modified  turbulence  model  were  determined  via 
numerical  optimization  for  benchmark  jet  flows  unrelated  to  the  impinging  jet  flows  of 
primary  concern. 

In  three  cases  considered,  two  impinging  jets  issuing  from  a  nozzle  normal  to  the 
ground  plane  and  one  issuing  from  a  nozzle  inclined  45°  into  the  crossflow,  the  modified 
model  consistently  improved  the  accuracy  for  the  computed  location  of  the  ground  vortex 
Errors  in  the  computed  vortex  location,  relative  to  the  distance  between  the  impingement 
point  and  the  ground  vortex  core,  were  as  large  as  50%  in  calculations  with  the  standard 
k— e  model  and  could  be  of  either  sign  (overprediction  or  underprediction  of  the  ground 
vortex's  forward  extent).  The  modified  model  reduced  these  errors  to  the  range  of  13%- 
18%,  and  it  consistently  overpredicted  the  upstream  penetration  of  the  ground  vortex. 
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1.  INTRODUCTION 


Because  of  recently  defined  tactical  threats  it  is  desirable  to  improve  the 
maneuverability  and  take-off  and  landing  performance  of  aircraft.  The  use  of  exhaust  jet 
thrust  vectoring,  including  thrust  reversing,  is  one  means  of  achieving  these  improvements. 
However,  thrust  vectoring  can  increase  significantly  the  complexity  of  the  flow  field  around 
the  aircraft,  and  make  aircraft  design  more  difficult.  To  achieve  good  aircraft  designs  in 
spite  of  this  difficulty,  it  is  necessary  to  have  an  accurate  flow  prediction  method.  The  flight 
regime  which  is  perhaps  the  most  difficult  to  predict  accurately  is  that  of  an  aircraft  landing 
with  thrust  reversers  operating.  The  flow  from  the  reversed  jets  creates  a  complex  flow 
environment  which  includes  an  upwash  and  ground  vortex,  for  example,  as  sketched  in 
Figure  1. 


Figure  1 .  Sketch  of  jet  with  thrust  reversers  in  ground  effect. 


The  effects  of  operating  thrust  reversers  in  ground  effect  have  been  investigated 
experimentally  in,  for  example.  References  1-3.  There  are  at  least  two  primary  problems 
caused  by  the  use  of  thrust  reversers.  One  is  that  forces  and  moments  acting  on  the  aircraft 
are  affected  significantly  by  the  reverser  plume.  The  other  is  that  hot  exhaust  gases  and 
debris  from  the  ground  may  be  thrown  up  near  the  front  of  the  aircraft  and  ingested  by  the 
engine,  which  causes  power  loss  and  engine  damage.  Both  the  ground  vortex  and  upwash 
can  be  highly  unsteady,  although  that  unsteadiness  is  usually  at  a  frequency  which  exceeds 
the  response  frequency  of  the  aircraft.  There  may  be  otner  problems  associated  with  thrust 
reversing,  for  example,  the  thermal  and  acoustic  loads  which  affect  the  aircraft  or 
undercarriage  stores. 

The  stability  and  control  problems  caused  by  thrust  reversing  depend  strongly  on  the 
aircraft's  proximity  to  the  ground.  Out  of  ground  effect,  20%-45 %  reductions  in  the 
effectiveness  of  the  horizontal  tail  have  been  noted,  while  the  effectiveness  of  the  vertical 
fins  and  ailerons  were  essentially  unaffected.  In  ground  effect  all  control  surfaces  are 
affected  1  Changes  in  aircraft's  lift  coefficient  of  the  order  of  +0.15/-0.9  and  changes  in 
pitching  moment  of  the  order  of  +0.5/-0.1  have  been  observed  by  Blake.3  These  forces 
depend  on  the  forward  deflection  and  outward  cant  of  the  reversing  jets,  the  ratio  of  jet  to 
freestream  dynamic  pressures,  on  the  height  above  the  ground,  and  the  side  slip  angle. 
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The  problem  of  hot  gas  and  debris  ingestion  has  long  been  a  research  issue  for  VTOL 
aircraft,4-5  because  engines  must  be  operated  near  maximum  power  for  the  aircraft  to  hover 
and,  hence,  in  danger  of  compressor  stall  if  hot  gases  are  ingested.  Both  the  ground  vortex 
and  the  upwash  fountain  may  contribute  to  contaminant  ingestion  by  causing  jet  exhaust 
gases  to  flow  up  and  in  front  of  the  aircraft.  On  an  aircraft  landing  with  thrust  reversers,  the 
primary  factors  which  determine  if  contaminant  ingestion  will  occur  are  the  pattern  and  the 
forward  extent  of  the  ground  vortex. 

Prediction  Methods  -  The  use  of  an  accurate  flow  prediction  method  based  on  the 
Navier-Stokes  equations  is  a  valuable  means  of  understanding  the  above  problems  and 
designing  aircraft  to  avoid  them.  Navier-Stokes  prediction  methods  are  being  used  in  a 
wide  variety  of  research  and  design  efforts  because  they  are  cost-effective,  fast  to  use  and 
can  provide  flow  details  which  are  often  not  obtained  in  experimental  studies.  Many 
experimental  programs  now  include  a  complementary  computational  program  as  a  means 
of  enhancing  the  value  of  the  experiments.  With  the  development  of  improved  computers 
and  numerical  algorithms,  the  accuracy  and  capabilities  of  these  prediction  methods  has 
increased  significantly  in  the  past  few  years.  Further  improvements  are  expected 
However,  there  is  one  critical  area  in  which  improvements  are  not  being  made  rapidly,  and 
that  is  turbulence  modeling. 

The  turbulence  generated  in  ground  effects  flows  is  far  more  complex  than  in 
aerodynamic  flow's  on  conventional  aircraft  at  moderate  incidence,  and  it  has  a  strong  effect 
on  the  complete  flow  field.  The  ground  vortex  and  upwash  fountain  are  strongly 
dependent  on  turbulent  mixing,  and  they  can  exert  a  strong  influence  on  aerodynamic 
characteristics  in  ground  effect. 

Turbulence  Modeling  -  Many  aspects  of  turbulence  in  ground  effects  flows  are  not 
predicted  well  by  existing  Reynolds-averaged  turbulence  models.  For  example,  the  size  and 
location  of  the  ground  vortex  is  predicted  poorly  by  the  Baldwin-Lomax  model  and  by  the 
standard  k— e  model.  The  shear  stresses  in  the  upwash  fountain  are  underpredicted  by  as 
much  as  an  order  of  magnitude  by  the  k-e  model,  the  Algebraic  Stress  Model  (ASM),  and 
the  Reynolds-Stress  transport  model  (RSTM).7  These  flows  are  highly  dependent  upon 
turbulent  momentum  transport,  and  the  errors  in  the  predicted  turbulence  create  large 
errors  in  the  mean  flow  predictions.  The  magnitude  of  these  modeling  errors  may  render 
existing  flow  prediction  methods  almost  useless  for  studies  of  the  thrust  reverser  effects  on 
landing  aircraft. 

The  prospects  for  developing  improved  turbulence  models  by  conventional  means  is 
not  good.  The  turbulence  mechanisms  which  are  important  in  these  flows  are  more 
complex  than  in  typical  aerodynamic  flows.  They  are  also  more  powerful  in  impinging 
flows,  meaning  that  they  produce  larger  and  more  rapid  changes  in  turbulence  than  in 
typical  aerodynamic  flows.  A  "turbulence  mechanism"  denotes  anything  that  affects  the 
turbulence  in  a  flow.  For  example,  the  conversion  of  mean  flow  kinetic  energy  into 
turbulence  kinetic  energy  due  to  the  interaction  between  turbulent  stresses  and  mean  flow 
gradients  is  one  such  mechanism.  Many  of  these  mechanisms  are  difficult  or  impossible  to 
identify  in  experimental  data,  and  they  can  have  a  strong  impact  on  these  flows.  Because  of 
the  rapid  spatial  variations  in  the  mean  flow,  turbulence  mechanisms  can  also  change 
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rapidly.  Some  modeling  concepts  which  are  widely  employed,  such  as  isotropy  and  the  use 
of  an  eddy  viscosity,  may  be  invalid  for  these  flows. 

The  technology  of  large  eddy  simulation  (LES)  provides  a  means  of  obtaining  the 
turbulence  data  which  are  required  to  identify  and  model  the  critical  turbulence 
mechanisms  in  these  flows.  LES  involves  numerical  simulations  of  fluid  flow,  based  on  the 
unsteady  Navier-Stokes  equations,  which  resolve  the  mean  flow  and  the  large  scales  of 
turbulence  Hence,  LES  can  represent  with  good  accuracy  the  dominant  features  of  a 
turbulent  flow.  All  turbulence  information  in  an  LES  can  be  recorded,  and  virtually  any 
turbulence  characteristics  can  be  studied.  The  detailed  knowledge  of  turbulence  provided 
by  LES  is  essential  to  identifying  and  modeling  the  complex  mechanisms  which  occur  in 
thrust  reverser  flows. 

Previous  studies  of  thrust  vectoring  in  ground  effect  have  focused  largely  on  vertical 
jets  used  with  Vertical  Take-off  and  Landing  (VTOL)  aircraft,  and  the  specific  problems  of 
thrust  reversing  have  not  been  addressed  on  a  fundamental  basis.  Therefore,  the  specific 
fluid  mechanical  problems  associated  with  thrust  reversing  are  not  well  understood.  Given 
the  limited  understanding  of  turbulence  in  impinging  jets,  it  is  reasonable  to  expect  that 
changing  the  inclination  of  the  jet  roughly  45®  may  have  a  significant  effect  on  the 
turbulerce.  Hence,  measurements  of  turbulence  in  thrust  reverser-like  flows  are  needed. 
One  study  which  included  inclined  jets8  involved  only  flow  visualization,  which  does  not 
provide  the  information  required  to  develop  improved  turbulence  models. 

Modeling  Priorities  -  In  a  previous  study  by  some  of  the  same  authors,  various 
aspects  of  the  types  of  flows  which  could  be  generated  in  thrust  reverser  applications  were 
evaluated.  The  following  list  assigns  priorities  to  the  generic  flow  phenomena  for  which 
improved  turbulence  models  are  most  critically  needed,  based  on  that  flow's  importance  to 
aircraft  safety  and  performance,  and  on  the  shortcomings  of  existing  turbulence  models. 

1.  The  ground  vortex  exerts  a  strong  influence  on  most  aerodynamic  characteristics, 
such  as  lift,  pitching  moment,  and  rolling  moment,  and  it  can  cause  contaminant 
ingestion.  The  position  of  the  ground  vortex  is  not  predicted  well  by  algebraic 
models  or  the  standard  k-£  model.  The  ground  vortex  flow  consists  of  the  vortex 
itself  and  other  closely  associated  regions,  the  wall  jet  and  the  impingement  zone. 

The  wall  jet  is  one  distinct  flow  region  of  the  ground  vortex.  The  forward  extent  of 
the  ground  vortex  is  largely  determined  by  the  deceleration  and  separation  of  the  wall 
jet.  The  wall  jet  beneath  the  ground  vortex  is  significantly  different  from  that  in  many 
laboratory  experiments  in  that  it  has  a  complex  vortical  flow  above  it,  and  it  is  highly 
three-dimensional.  The  accuracy  with  which  this  region  of  flow  can  be  computed  is 
not  known. 

The  impingement  zone  is  where  the  jet  strikes  the  ground  and  spreads  out  into  a  wall 
jet.  The  jet  flow  must  pass  through  this  region  before  reaching  the  upwash  and  the 
ground  vortex,  so  that  inaccuracy  here  will  propagate  to  most  of  the  ground  effects 
region.  This  region  is  predicted  poorly  by  the  standard  k-c  model,  although  it 
appears  that  means  for  significantly  improving  the  accuracy  in  this  region  exist. 
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2.  The  upwash  fountain  and  related  flow  phenomena  tend  to  throw  contaminants  up 
and  in  front  of  an  aircraft  where  they  can  be  ingested.  In  experiments  this  problem 
does  not  appear  to  be  as  serious  as  the  ground  vortex.  However,  various  Navier- 
Stokes  calculations,  performed  at  Nielsen  Engineering  &  Research  (NEAR), 
NASA/ Ames  Research  Center,  and  NASA/ Lewis  Research  Center,  all  predict  the 
upwash-like  phenomena  to  be  more  significant  than  the  ground  vortex  for  causing 
ingestion  problems.  This  result  is  not  supported  by  experiment.  Accurate  modeling 
of  the  upwash  and  related  phenomena  is  required  to  reduce  its  predicted  effect  to  a 
realistic  magnitude. 

3.  A  free  jet  in  a  crossflow  is  the  first  phenomenon  to  occur  as  the  flow  emerges  from  the 
thrust  vectoring  apparatus.  As  the  crossflow  interacts  with  the  jet,  two  counter 
rotating  vortices  are  formed  on  the  sides  of  the  jet.  If  the  jet  is  adjacent  to  a  vertical 
stabilizer,  the  vortices  combined  with  the  blockage  effect  of  the  jet  can  significantly 
alter  the  pressure  acting  on  the  fin.  Turbulence  in  vortices  is  not  modeled  well  by 
existing  models.  The  effects  of  the  mean  and  turbulence  quantities  at  the  jet  exit  can 
play  a  significant  role  in  the  jet's  development. 

4.  Dynamical  effects,  due  to  the  rate  of  descent  of  the  aircraft,  are  important.  It  is  not 
currently  known  if  modeling  of  dynamical  effects  on  turbulence  is  a  research  issue. 
The  time  scales  for  turbulence  are  much  shorter  than  for  the  mean  flow,  so  that  a 
quasi-steady  model  may  suffice.  Also,  the  accuracy  for  steady  state  flow's  is  very 
poor,  at  present,  and  must  be  improved  before  an  attempt  at  modeling  unsteady  flow 
is  appropriate. 

The  ground  vortex  and  upwash  are  the  most  critical  to  be  able  to  predict  accurately.  If 
they  can  be  predicted  accurately,  it  is  likely  that  the  other  phenomena  can  be  predicted 
because  they  are  governed  by  many  of  the  same  turbulence  mechanisms  which  occur  in  the 
first  two. 

Objectives 

The  objectives  of  this  work  are  to  develop  and  to  validate  a  Reynolds-averaged 
turbulence  model  for  thrust  reverser  flow  fields.  Successful  completion  of  this  work  will 
yield  an  accurate  mathematical  model  of  turbulence  which  can  be  incorporated  into  Navier- 
Stokes  prediction  codes  for  application  to  thrust  reverser  problems.  To  achieve  these  goals, 
several  intermediate  objectives  must  be  met. 

The  first  objective  is  to  obtain  good  turbulence  data  for  these  flows  so  that  the 
mechanisms  which  control  turbulence  can  be  evaluated.  Data  will  be  obtained  from 
experiments  and  from  large  eddy  simulations  (LES)  of  impinging  jet  flows.  Experiments  are 
needed  because  they  provide  trustworthy  data  which  is  accurate  to  within  some 
experimental  uncertainty.  However,  some  important  turbulence  data,  for  example,  two- 
point  correlations  and  fluctuating  pressures,  are  difficult  or  impossible  to  obtain  from 
physical  experiments,  at  present.  Therefore,  LES  will  also  be  performed  to  obtain  data 
which  are  currently  inaccessible  in  experiments. 

The  next  objective  is  to  identify  the  mechanisms  which  control  turbulence  in  the 
various  flow  regions.  Streamline  curvature,  anisotropy,  and  large  scale  mixing  are 
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important  mechanisms  in  the  upwash,  and  they  are  expected  to  be  important  in  other 
regions,  as  well.  These  are  fundamental  mechanisms  which  are  somewhat  understood,  and 
which  can  be  identified  in  the  LES  results.  There  are  also  other  mechanisms  which  are  not 
yet  well  understood.  For  example,  ihe  position  of  the  ground  vortex  is  determined  by  the 
deceleration  and  separation  of  a  wall  jet  beneath  a  vortical  outer  flow  moving  in  the 
opposite  direction.  The  mechanisms  which  control  this  process  have  not  yet  been  identified, 
but  they  must  be  understood  if  a  model  which  can  accurately  predict  the  forward 
penetration  of  the  ground  vortex  is  to  be  developed. 

The  third  objective  is  to  develop  modifications  to  the  k— e  model  which  account  for  the 
dominant  additional  turbulence  mechanisms  associated  with  jet  impingement.  Based  on  the 
understanding  of  turbulence  mechanisms,  a  model  or  models  will  be  developed  which 
mimic  the  turbulence  behavior  observed  in  the  experiments  and  LES.  Model  development 
is  one  of  the  most  critical  steps  in  the  proposed  work.  Many  individual  mechanisms  are 
poorly  understood  at  present  and  predicted  poorly  or  not  at  all  by  most  turbulence  models. 
The  model  must  also  be  optimized  to  give  the  best  possible  accuracy  for  the  range  of  flow 
phenomena  which  occurs  in  thrust  reverser  flow  fields.  This  will  involve  "fine  tuning"  the 
model  to  give  good  agreement  with  LES  and  experimental  data. 

The  final  objective  is  to  validate  the  accuracy  of  the  model  which  is  developed.  The 
turbulence  model  will  have  been  developed  and  optimized  for  the  specific  phenomena 
which  occur  in  a  few  critical  regions  of  thrust  reverser  flows.  It  is  important  to  determine 
how  well  the  model  performs  for  complete  experimental  flow  fields.  If  the  model  does  have 
significant  shortcomings  for  the  complete  flow  field,  it  will  be  necessary  to  tune  the  model 
further  to  give  the  best  possible  overall  accuracy.  If  significant  errors  still  persist,  it  is 
essential  to  know  where  they  occur  and  how  large  they  are. 
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2.  EXPERIMENTAL  STUDY 

This  part  of  the  project  was  performed  in  the  Mechanical  Engineering  Department  of 
Stanford  University  by  Prof.  P.  Bradshaw  and  several  research  assistants. 

Although  the  target  flow  is  a  hot,  transonic  jet  moving  relative  to  the  ground,  the  main 
problems  in  turbulence  modeling  will  occur  in  a  low-speed,  isothermal  flow  in  which  the  jet 
is  stationary  relative  to  the  ground.  The  results  will  be  quantitatively  different  from  the  hot, 
high-speed  flows  resulting  from  engine  exhaust  jets.  However,  data  from  the  low-speed 
case  can  provide  a  good  understanding  of  the  qualitative  characteristics  of  thrust  reverser 
flows. 

Wind  Tunnel  -  The  tunnel  used  for  this  experiment  is  an  open-circuit,  suck-down 
tunnel  designed  especially  for  smoke-flow  visualization.  The  intake  consists  of  a  30  in  x  30 
in  bell-mouth  followed  directly  by  a  flow-management  section.  In  this  section,  the  incoming 
flow  passes  through  a  1  in  long  honeycomb  straightener  of  0.5  in  cell  diameter  and  then 
through  a  set  of  four  screens.  The  screens  are  separated  by  1.5  in  and  each  has  a  pitch  of  38 
wires? inch.  The  wire  diameter  is  0.006  in  which  translates  to  an  open  area  ratio  of  0.5%. 
There  is  no  contraction  to  the  working  section  which  begins  immediately  downstream  of  the 
screens.  A  diagram  of  the  tunnel  is  given  in  Figure  2. 


The  working  section  floor  is  made  from  particle  board  which  has  been  painted  black  to 
reduce  glare  during  flow  visualization.  The  walls  are  made  from  1  /4  in  Plexiglass,  and  the 
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roof  is  formed  by  interchangeable  particle-board  panels.  Although  constant  cross  section  is 
maintained  for  a  length  of  120  in,  the  end  of  the  working  section  is  effectively  defined  by  the 
region  where  interference  from  the  exit  begins  (roughly  90  in  from  the  start  of  the  test 
section). 

The  pressure  drop  through  the  tunnel  is  supplied  by  a  20.8  in  intake-diameter, 
centrifugal  blower  driven  by  a  constant  speed  motor.  For  this  experiment,  the  maximum 
flow  rate  of  16  ft/s  is  desired.  The  flow  rate  is  measured  by  a  pitot-static  port  downstream 
of  the  flow-management  section. 

Cross-flow  Jet  -  An  auxiliary  blower  is  used  to  power  the  cross-flow  jet.  The  intake  to 
this  blower  is  restricted  by  a  perforated  aluminum  sheet  and  a  butterfly  valve  which  can  be 
adjusted  to  maintain  a  constant  ratio  of  bulk-average  jet  speed  to  tunnel  speed  of  4.36.  The 
outflow  is  led  through  a  4.25  in  diameter  flexible  duct  to  a  contraction  which  leads  into  a  2 
in  inside  diameter  feed  pipe.  At  the  entrance  to  this  pipe  is  a  2.5  in  section  of  aluminum 
honeycomb  with  a  1/4  in  cell  diameter.  The  feed  pipe  passes  through  the  roof  at  45*  and 
extends  23  in  into  the  tunnel.  The  end  of  this  pipe  is  supported  by  an  aluminum  "sting" 
which  is  supported  by  aluminum  rods  fastened  to  the  floor.  A  conical  contraction  takes  the 
flow  into  the  1.131  in  inside  diameter  nozzle  section. 

The  nozzle  section  extends  for  16  in  ending  with  its  exit  plane  inclined  at  45*  and  its 
center  3.5  in  above  the  working  section  floor  and  39  in  downstream  of  the  final  screen.  (See 
Figure  2.)  While  a  contraction  at  the  nozzle  exit  could  have  provided  uniform  flow,  it  was 
felt  more  important  to  minimize  the  nozzle  area  and  thus  its  interference  with  the  flow. 
Furthermore,  the  flow  exiting  actual  thrust  reversers  will  probably  be  highly  turbulent 
because  of  the  engine  turbulence  and  the  turning  vanes.  Thus,  the  reversed  engine  exhaust 
is  more  closely  approximated  by  a  turbulent  pipe  flow  than  a  jet  from  a  large  area-ratio 
contraction.  The  flow  rate  is  measured  by  a  pitot  tube  in  the  contraction  to  the  feed  pipe 
and  a  static  port  3  in  from  the  exit  of  the  nozzle.  A  Plexiglass  template  was  manufactured  to 
check  the  alignment  of  the  jet  nozzle. 

Tunnel  Calibration  -  A  series  of  tests  were  run  to  determine  the  quality  of  the  flow 
entering  the  working  section  of  the  tunnel.  The  first  test  was  a  total  pressure  survey  3  in 
downstream  of  the  flow-management  section.  An  initial  traverse  indicated  significant  non¬ 
uniformity  in  the  velocity  field.  The  results  of  a  second  survey,  performed  after  cleaning  the 
screens,  are  shown  in  Figure  3,  and  represent  a  marked  improvement.  In  order  to  keep  dirt 
and  dust  interference  from  altering  the  flow  field,  the  screens  must  be  cleaned  after  every  2 
hours  of  tunnel  operation. 
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Figure  3.  Velocity  uniformity  at  inlet  to  wind  tunnel  test  section  from  total  pressure 
survey. 


The  second  diagnostic  test  was  a  series  of  total  pressure  traverses  of  the  floor  boundary 
layer.  These  were  done  on  the  tunnel  center-plane  at  distances  of  7  and  16  inches  from  the 
final  screen,  with  a  1  /8  in  high  trip  of  square  cross  section  placed  3  in  downstream  of  the 
screens.  The  first  profiles  revealed  a  surprisingly  large  boundary  layer  thickness.  An 
investigation  was  done  to  determine  if  this  was  a  result  of  dirty  screens,  intake  geometry,  or 
trip  shape.  Comparing  a  traverse  done  without  the  trip  to  one  made  of  the  roof  boundary 
layer  suggested  the  trip  was  the  culprit.  Replacing  the  square  trip  with  a  cylinder  of  equal 
height,  faired  by  a  tape  strip,  resulted  in  the  boundary  layer  shown  in  Figure  4.  This  layer  is 
felt  to  be  of  acceptable  thickness.  The  "wiggles"  near  the  top  of  the  layer  do  not  appear  in 
traverses  done  without  the  trip  and  are  felt  to  be  a  consequence  of  the  trip. 
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Jet  Calibration  -  Vertical  and  horizontal  traverses  of  the  jet  were  performed  at  a 
station  0.55  in  inside  the  nozzle.  The  results  of  these  traverses  are  shown  in  Figure  5-a  along 
with  computed  mass  and  momentum  flow  rates.  In  order  to  compare  with  the  "law  of  the 
wall"  profile,  shear  stress  measurements  were  made  at  the  nozzle  wall  and  a  logarithmic 
plot  of  u+  vs.  y+  was  generated.  This  plot,  along  with  the  law  of  the  wall  profile,  is  shown  in 
Figure  5-b. 


Ja  Dyqwik  Prcnure  Profit 


(a)  Nozzle  exit  velocity  profile 
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J«  Velocity  Profile  (tof  ptot) 


(b)  Near-wall  velocity  profile 
Figure  5. 


Selection  of  the  Primary  Test  Case  -  In  order  to  generate  a  flow  field  representative  of 
a  thrust  reverser,  a  series  of  flow  visualization  experiments  were  done  to  select  the  angle  of 
the  nozzle  and  its  height  above  the  floor.  With  the  flow  into  the  auxiliary  blower  fully 
restricted,  the  nozzle  rotated  to  0  =  43.3°  forward  of  vertical  and  the  center  of  the  nozzle 
3.57  in  above  the  floor,  the  flow  pattern  shown  in  Figure  6  resulted.  The  final  geometry  was 
fixed  very  near  this  setting  because  it  was  decided  that  this  configuration  provided  a 
representative  recirculation  of  the  reversed  flow.  The  actual  geometry  shown  in  Figure  2 
was  selected  because  it  provides  round  figures  for  the  test  case. 


Figure  6.  Smoke  visualization  photograph  of  inclined  impinging  jet. 
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The  primary  inclined  jet  test  used  in  this  experiment  and  computed  in  Reynolds- 
averaged  and  large  eddy  simulation,  described  later,  is  defined  by  the  following  operating 
conditions: 


0  =  45° 

U_  =  16  ft. /sec 
Hj/Dj  =  3.09 

Vj,tnax/Vj,bulk=  1208 

Vtmax/U.  =  5.266 


nozzle  angle  forward  of  straight  down 
crossflow  velocity 
nozzle  height  to  diameter  ratio 
maximum  to  bulk  jet  velocity  ratio 
jet  maximum  to  crossflow  ratio 


Boundary  Conditions  -  In  order  to  model  the  flow  computationally,  appropriate 
boundaries  must  be  selected  and  the  flow  field  determined  there.  For  this  flow,  the  in-flow 
boundaries  are  a  tunnel  cross  section  downstream  of  the  trip  and  the  nozzle  exit  plane. 
These  measurements  should  be  performed  with  both  the  tunnel  and  the  jet  running. 
However,  the  tunnel  inflow  conditions  were  measured  with  the  jet  turned  off.  These 
include  the  boundary  layer  traverses  already  mentioned.  The  jet  exit  conditions  were 
measured  with  a  three  tube  yawmeter  with  the  tunnel  running.  The  purpose  of  the  jet 
traverse  was  to  determine  how  much  the  free-steam  deflects  the  jet  exit  flow.  The  results  of 
this  traverse  are  shown  in  Figure  7. 


Figure  7.  Deflection  angles  at  jet  exit. 

Floor  Static  Pressure  Distribution  -  The  static  pressure  distribution  along  the  center- 
line  of  the  floor  was  measured  by  means  of  an  array  of  pressure  taps  built  into  the  floor. 
The  spacing  of  these  taps  was  chosen  to  provide  a  high  density  of  measurements  in  the 
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regions  of  rapidly  changing  pressure.  The  results  of  the  static  pressure  traverse  for  the 
primary  test  case  are  shown  in  Figure  8. 


Floor  Static  Pressure  Distribution 


Figure  8.  Surface  static  pressure  distributions  for  primary  test  case. 


To  provide  data  for  model  testing,  the  jet  velocity  was  varied  over  a  wide  range  while 
the  jet  height  and  angle  were  held  constant.  The  surface  static  pressures  are  given  in 
Figure  9  for  several  jet  velocities.  The  ground  vortex  does  move  forward  with  increasing  jet 
velocity,  and  there  are  changes  in  the  shape  of  the  pressure  distribution  as  the  velocity  ratio 
changes.  Note  that  the  scale  of  these  figures  changes  to  accommodate  the  changing  velocity 
ratio. 
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Figure  9.  Surface  static  pressure  distributions  for  several  velocity  ratios. 

Floor  Shear  Stress  Distribution  -  The  shear  stress  distribution  was  measured  in  the  same 
region  as  the  static  pressure  distribution  by  means  of  a  Preston  tube  traverse.  The  outside 
diameter  of  the  Preston  tube  used  was  0.127  in  and  the  calibration  used  was  that  of  Patel.9 
In  order  to  make  shear  stress  measurements  in  the  reversed  flow  regions,  a  complete 
traverse  was  made  with  the  probe  facing  the  tunnel  entrance,  followed  by  a  similar  traverse 
with  the  probe  direction  reversed.  The  pressure  difference  at  any  given  location  was  taken 
to  be  the  greater  of  the  two  measurements.  The  shear  stress  distribution  appears  in 
Figure  10.  Note,  however,  that  this  data  corresponds  to  a  flow  configuration  of  $  =  43.3° 
and  H  /D  =  3.16,  which  is  slightly  different  from  the  primary  test  case. 
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Figure  10.  Surface  centerline  skin  friction. 
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3.  NUMERICAL  METHODS 

The  focus  of  this  study  is  turbulence  rnd  its  modeling.  However,  much  of  the  work  has 
been  conducted  using  numerical  solutions  of  the  governing  flow  equations,  and  the  relevant 
features  of  the  numerical  simulation  methodology  are  given  here. 

Navier-Stokes  Equations 

All  calculations  described  in  this  report  are  based  upon  numerical  solutions  of  the 
Navier-Stokes  equations  in  finite-volume  form.  The  discretized  solution  is  represented  by 
Q,  the  average  value  in  a  cell  volume  V,  defined  by 

Q  =  -  [  Q  dV  (1) 

JV 

in  which 

Q  =  (p,pu,pv,p  w,pe)‘  (2) 

is  the  vector  of  conserved  variables,  mass,  the  momenta,  and  total  energy.  The  conservation 
law  equations  for  these  variables,  the  "Navier-Stokes"  equations,  can  be  written  as 


[  P-dS  =  0 
JS 


(3) 


in  which  the  surface  area  of  the  finite  volume  is  denoted  S,  dS  is  the  outward  surface 
normal,  and  F  is  the  flux  tensor  of  the  conserved  variables: 
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The  thermodynamic  pressure  is  given  by 


p  =  (y-1) (pe  -  0.5p(u2  +  v2  +  w2)) 


(5) 


in  which  e  =  cvT  +  (u2+v2+w2)/2  is  the  total  internal  energy.  The  viscous  and  turbulent 
stresses  and  the  thermal  stresses  are  given  in  Cartesian  tensor  notation  by 


duA  du 

Tij  -  <"  +  "t)(  +  sq 


(k  *  V  H 


[6) 


(7] 


j 


in  which,  for  example,  t12  =  tx  ,  T  is  temperature,  nt  is  a  modeled  turbulent  viscosity,  and  k, 
is  a  modeled  turbulent  thermal  diffusivity. 

Equations  (1-7)  which  describe  the  fluid  motion  are  exact,  even  when  applied  on  a  grid 
with  discrete  finite-volume  cells.  The  interpolation  and  differentiation  schemes  used  to 
compute  the  fluxes  on  the  cell  faces  determine  the  accuracy  of  the  spatial  discretization 
scheme.  These  are  described  below. 

Spatial  Discretization  Scheme 

High-accuracy,  central,  discretization  algorithms  are  described  for  the  finite-volume 
form  of  the  Navier-Stokes  equations  in  this  section.  The  higher  order  forms  of  these 
methods  were  developed  recently,10 and  no  description  of  them  has  yet  been  published. 
Hence,  a  summary  of  the  methods  is  given  here. 

The  integral  or  finite  volume  form  of  the  Navier-Stokes  equations,  Eq.  (3)  is  exact; 
however,  it  becomes  approximate  when  numerical  approximations  are  used  to  compute  the 
fluxes  through  the  faces  of  the  finite-volume  cells.  Interpolation  is  used  to  compute  the 
inviscid  fluxes,  while  finite  difference  approximations  are  used  to  compute  viscous  fluxes. 
The  inviscid  fluxes  are  the  more  important  of  these  two,  for  reasons  discussed  below. 

High-Accuracy  Spatial  Interpolation  -  The  accuracy  of  the  spatial  scheme  depends  on 
the  accuracy  of  the  interpolation  used  to  determine  the  values  of  the^cell-face  fluxes,  based 
on  the  volume-average  value  of  the  conserved  variables  in  the  cells.  Q  is  not  associated  with 
any  specific  location  in  the  cell,  and  specifically,  Q  is  not  generally  the  value  of  Q  at  the  cell 
center.  This  fact  must  be  recognized  when  developing  the  discretization  scneme. 
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Co-location  interpolants  which  pass  through  Q  assigned  to  the  cell  centers  fail  to  heed  this 
warning,  and  they  are  at  best  second  order  accurate.  An  interpolation  scheme  based  on  the 
assumption  that  the  value  at  the  cell  center  is  the  average  value  will  incur  an  error 
proportional  to  Ax2d2F/dx2.  This  is  the  same  magnitude  error  as  incurred  by  a  linear 
interpolation  and  does  not  degrade  the  error  of  a  second  order  algorithm.  If  a  higher  order 
interpolation  scheme  is  used,  it  is  necessary  to  account  for  the  nonlinear  distribution  of  the 
solution  within  the  cell. 


In  the  following  equations,  a  fourth  order  accurate  interpolation  scheme  will  be 
developed  in  one_space  dimension.  Fourth  order  is  the  lowest  order  for  which  special 
recognition  that  Q  does  not  coincide  with  the  cell  center  is  required.  Following  that, 
interpolation  formulae  will  be  given  for  schemes  up  to  eighth  order. 

The  Navier-Stokes  equations  can  be  written  in  one  space  dimension  as 


at  Q  +  V  [Fj+l/2“  Fj-l/2]  = 


(8) 


in  which  A  is  the  area  of  the  faces  of  the  finite  volume,  Fj+1/2  and  Fj_]/2  are  the  fluxes  at  the 
cell  faces,  and  V  is  the  cell  volume.  Equation  (8)  is  exact,  as  is  Eq.  (3);  approximation  is 
introduced,  and  the  spatial  accuracy  of  the  scheme  is  determined  by  the  method  used  to 
obtain  the  cell  face  fluxes  Fi+)/j  from  Q. 


A  simple  means  of  achieving  this  accuracy  is  to  require  that  the  difference  of  the 
interpolants  across  a  cell  be  equal  to  a  finite  difference  approximation  of  the  appropriate 
order.  Hence, 
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(9) 


must  possess  desired  accuracy,  where  fiF/  fix  is  the  finite  difference  approximation  with  the 
desired  order  of  accuracy.  From  the  fourth  order  finite  difference  approximation, 

fiF /fix  =  (F^  -  8FH  +  8Fj+1  -  Fj+2)/12Ax  (10) 

the  fourth  order  interpolation  is  found  to  be 

Fj+i,j  =  +  7Fj  +  7Fj+i  *  Fj+2>/12  (11) 

In  contrast,  the  fourth  order  accurate  co-location  interpolant  which  wrongly  assumes  that  Q 
exists  at  the  cell  center  is  given  by 

Fj+t,j  -  £'Fj-i  +  9Fj  +  9Fj+,  -  Fj+2)/16  (12) 


When  differenced,  this  interpolant  yields  the  approximation  of  the  gradient 

fiF/ fix  =  (F|  2  -  10FH  +  10Fj+1  -  Fj+2)/l6Ax  (13) 
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whose  leading  error  term  is  0(Ax2).  Hence,  Eq.  (11)  gives  the  true  fourth  order  accurate 
representation  of  the  flux  at  a  cell  face,  which  accounts  for  the  variation  of  Q  within  the 
cell's  volume. 


Formulae  for  high-order  interpolation  can  be  developed,  following  the  approach  in 
Eqs.  (9-11),  using  Taylor  series  expansions.  Central  interpolation  schemes  through  eighth 
order  are  given  in  Table  1 . 


Table  1.  Interpolation  Weights  for  Fj+J/ 
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These  formulae  are  based  on  the  assumption  of  uniform  cell  sizes,  and  expressions  for 
stretched  grids  have  not  yet  been  developed,  ^rid  stretching  does  reduce  the  accuracy  of 
the  interpolation  schemes.  Therefore,  the  grid,  used  in  the  present  study  are  as  smooth  as 
practically  possible  in  regions  where  the  solution  must  be  computed  with  good  accuracy. 

Spatial  Accuracy  in  Fourier  Space  -  The  principal  reason  for  using  a  discretization 
method  with  a  high  order  of  accuracy  is  to  minimize  the  number  of  grid  points  needed  to 
accurately  represent  variations  in  the  solution.  Hence,  a  very  relevant  means  of  assessing 
the  accuracy  of  a  method  is  to  determine  how  few  grid  points  are  required  to  accurately 
characterize  some  aspect  of  the  solution.  Fourier  analysis  provides  a  good  means  of 
performing  this  evaluation. 

The  shortest  wave  length  represented  on  a  grid  with  a  spacing  Ax  is 

/l0  =  2Ax  (14) 

which  corresponds  to  oscillations  at  every  other  grid  point.  The  equivalent  highest  wave 
number  is 


(Jq  =  rr/Ax  (15) 

A  numerical  discretization  scheme  will  give  accurate  derivatives  up  to  some  "cutoff’ 
wave  number  uc  <  <a0,  beyond  which  the  accuracy  falls  rapidly.  wc  depends  on  the 
derivative  being  computed  (e.g.,  first,  second),  the  level  of  accuracy  required,  and,  of  course, 
the  discretization  scheme. 
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Figure  11.  Accuracy  of  spatial  discretization,  defined  in  Eq.  16. 


The  results  of  a  Fourier  analysis  of  the  schemes  represented  in  Table  1  are  given  in 
Figure  11. 10  A  trial  function  F  =  sin (wx* /Ax)  is  used.  The  measure  of  accuracy  is: 


Accuracy  = 


6F/6x 

aF/dx 


(16) 


in  which  aF/dx  and  SF/ix  are  the  values  of  the  derivative  computed  exactly  and  with  the 
discretization  formulae.  Derivative  terms  are  computed  from  the  interpolation  formulae 
using  Eq.  (9).  All  interpolation  schemes  in  Table  1  have  100%  accuracy  as  w/ci^  ■*  0  and  zero 
accuracy  as  w/Wq  -»  1.  However,  the  higher  order  schemes  have  nearly  perfect  accuracy  to  a 
much  higher  wave  number  than  the  lower  order  schemes.  The  second  order  scheme  is 
better  than  99%  accurate  up  to  a  cutoff  wave  number  of  about  wc  =  0.05  (Jq,  while  the  eighth 
order  scheme  is  99%  accurate  to  about  uc  =  0.4  % 


In  this  comparison,  the  eighth  order  scheme  requires  roughly  1  /8th  as  many  grid 
points  in  each  spatial  direction  for  the  same  accuracy  as  the  second  order  scheme,  or  1  /500th 
as  many  grid  points  in  a  three-dimensional  calculation. 

Artificial  Dissipation  -  The  central  discretization  schemes  described  above  re  non- 
dissipative.  That  is,  energy  cannot  be  dissipated  by  numerical  discretization  of  the  inviscid 
fluxes.  These  schemes  permit  high  wave  number  oscillations  in  the  solution,  which  are  non¬ 
physical  and  can  cause  the  solution  to  diverge  if  they  are  not  regulated.  An  eighth  order 
nonlinear  dissipation  scheme  is  used  to  remove  energy  from  the  high  wave  numbers  in  the 
solution.  Because  the  dissipation  is  of  high-order,  low-wave  numbers,  roughly  <o  <  0^/2,  are 
unaffected  by  the  dissipation  and  computed  with  good  accuracy.  This  nonlinear  scheme 
was  developed  in  conjunction  with  the  high-order  discretization  scheme10  for  numerical 
problems  in  computational  aeroacoustics.  The  present  simulations  deal  with  subsonic 
flows,  and  the  shock-capturing  capabilities  of  the  dissipation  are  not  needed. 


Viscous  Fluxes  -  The  viscous  fluxes  at  the  cell  faces  are  computed  with  second  order 
central  difference  approximations.  The  specific  difference  approximation  used  depends 
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upon  the  flux  required  and  its  location  on  the  grid  cell.  Let  the  subscripts  j  and  k  denote  the 
indices  in  the  x-  and  y-directions.  Then,  for  example,  the  x-gradient  term  on  the  cell  face 
between  the  j  and  j+1  cells  is  computed  as 

6u/6x 1  j+I/j  k  =  (Uj+i  k  -  Ujjjl/Ax  (17) 

A  y-derivative  term  computed  at  the  same  location  is  more  complex  because  it  must  be 
computed  at  the  j  and  j+1  cell  centers  and  then  averaged  to  the  cell  face. 

5u/6y !  j+l/i  =  (u^,  -  u^.,  +  Uj+jjk+i  -uj+u.1)/2Ay  (18) 

Equation  (17)  is  more  accurate  than  (18)  because  the  difference  is  taken  over  a  shorter 
interval.  In  general,  the  most  important  viscous  fluxes  are  computed  with  Eq.  (17). 

The  accuracy  of  Eq.  (17),  when  substituted  into  the  Navier-Stokes  equations,  has  been 
evaluated  using  Fourier  analysis,  as  described  in  the  previous  section.  The  second  order 
viscous  terms  are  relatively  close  in  accuracy  to  the  fourth  order  approximation  of  the 
inviscid  fluxes  in  the  lower  half  of  the  wave  number  spectrum.  Hence,  the  second  order 
viscous  fluxes  are  comparable  in  accuracy  to  the  fourth  order  for  inviscid  fluxes. 

The  fourth  order  interpolation  scheme  was  used  for  the  Reynolds-averaged 
calculations,  because  comparable  accuracy  in  the  viscous  and  inviscid  terms  is  appropriate. 
The  eighth  order  scheme  was  used  for  the  large  eddy  simulations  which  are  convection- 
dominated. 

Boundary  Conditions  -  A  wide  variety  of  boundary  conditions  are  required  to 
perform  the  range  of  calculations  in  this  study.  In  general,  standard  boundary  condition 
methods  are  used,  and  a  relatively  brief  description  of  these  is  provided  here.  Although  the 
calculations  are  performed  withjhe  finite  volume  method,  the  boundary  conditions  are 
imposed  on  the  solution  variable  Q  rather  than  on  the  fluxes. 

The  impingement  wall  is  treated  with  no-slip  adiabatic  conditions,  and  therefore, 

u  =  v  =  w  =  dp/d  y  =  dp/dy  =  0  (19) 

in  which  y  is  normal  to  the  wall.  One-sided,  second  order  approximations  are  used  to 
compute  the  derivatives  at  the  wall,  so  the  boundary  values  are  computed  with  an 
extrapolation  expression, 

p,  =  1.333p2-0.333p3  (20) 

in  which  the  subscript  ()t  denotes  the  boundary  values  and  others  are  interior  values. 
However,  the  shear  stress  between  the  first  and  second  grid  points,  effectively  the  wall 
shear  stress,  is  computed  with  "wall  functions."  A  logarithmic  law-of-the-wall  velocity 
profile  is  assumed  to  exist  between  the  wall  and  the  first  interior  grid  point,  and  the  shear 
stress  is  computed  to  be  consistent  with  this  profile.  The  velocity  at  the  first  interior  point  is 
restrained  only  by  this  shear  stress,  and  not  by  the  artificial  dissipation;  hence,  that  velocity 
can  be  large.  Tunnel  side  and  top  walls  were  treated  with  an  impermeable  slip  boundary 
condition. 


The  non-reflecting  boundary  conditions  developed  by  Thompson,11  with  some 
modifications,  are  used  for  the  non-wall  boundaries,  which  included  the  inflow,  outflow, 
and  occasionally  the  top  boundaries.  These  are  briefly  summarized  here,  using  finite 
difference  notation  for  compactness.  The  equations  for  one-dimensional  flow  can  be  written 
in  nonconservative  form  as 


dQ  »  dQ 

+  A  ~  =  0 

dt  ax 


(21) 


in  which  A  is  the  Jacobian  matrix  of  the  inviscid  flux  vector.  Through  eigenvector- 
eigenvalue  manipulations,  the  A  matrix  can  be  split  into  two  matrices  which  have  positive 
and  negative  eigenvalues,  A  =  A+  +  A'.  They  represent  the  information  which  propagates  in 
the  +x  and  -x  directions,  respectively.  Upwind  differencing  consistent  with  the  direction  of 
signal  propagation  is  then  used  to  compute  spatial  gradients,  and  Eq.  (21)  becomes 
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Under  the  assumption  that  no  information  enters  the  domain  through  a  boundary,  it  is 
appropriate  to  discard  the  A'  term  at  the  xmax  boundary  and  the  A+  term  at  the  xmin 
boundary  because  these  terms  represent  incoming  information.  The  boundary  conditions 
are  applied  by  integrating  Eq.  (22)  with  either  A+  or  A'  set  to  zero  at  the  boundary  points 
plus  the  full  spatial  gradient  terms  in  the  other  spatial  directions. 

It  was  necessary  to  modify  the  Thompson  boundary  conditions  at  several  boundaries. 
With  no  incoming  information  it  is  not  possible  to  control  the  ambient  pressure.  This 
difficulty  is  resolved  by  letting  the  boundary  pressure  relax  slowly  to  the  selected  ambient 
pressure.  At  each  time  step,  the  boundary  pressure  was  typically  made  up  of  99%  of  the 
pressure  computed  with  Thompson's  boundary  conditions  and  1%  of  the  ambient  pressure 
at  each  application  of  the  boundary  conditions.  Hence,  a  characteristic  time  for  the 
relaxation  to  ambient  pressure  was  100  time  steps.  Since  this  change  in  pressure  is 
physically  caused  by  incoming  acoustical  waves  which  are  isentropic,  the  density  is  adjusted 
to  preserve  the  entropy  at  the  boundary.  The  velocities  at  a  subsonic  inflow  boundary  are 
also  imposed  by  computing  the  boundary  velocities  as  a  weighted  average  of  the  values 
computed  with  the  Thompson  boundary  conditions  and  the  imposed  mean  values,  with 
similar  weights  to  those  used  for  the  pressure.  Hence,  the  inflow  velocity  could  deviate 
from  the  desired  value  during  transients  in  the  calculation,  but  it  would  asymptote  to  the 
correct  value  in  the  steady-state  solution. 

Turbulent  Inflow  Boundaries  -  The  boundary  layer  flow  at  inflow  boundaries  was 
specified  from  experimental  data,  when  available.  Typically,  the  only  data  available  from 
the  freestream  boundary  layers  were  mean  velocity  profiles  or  total  thickness.  Data  from 
Cimbala12  indicate  that  the  dominant  characteristics  of  impinging  jet  flows  depend  only 
weakly  on  the  boundary  layer  thickness.  In  the  calculations,  a  y‘/7  velocity  profile  was  used 
for  the  mean  velocity,  and  the  turbulence  was  assumed  to  be  that  of  a  flat  plate  boundary 
layer.  The  turbulence  energy  k  was  taken  from  typical  flat  plate  data,  and  t  was  set  to  give 
the  correct  shear  stress,  given  k  and  the  velocity  profile. 


For  the  primary  test  case,  the  jet  flow  issued  from  a  long  constant-diameter  pipe  and 
was  nearly  fully  developed.  The  measured  velocity  profile,  including  the  measured 
deflection  from  straight  flow,  was  used.  Turbulence  was  assumed  to  be  from  a  fully 
developed  pipe  flow  at  the  appropriate  Reynolds  number. 

Large  Eddy  Simulation  Inflow  Conditions  -  Boundary  conditions  used  for  LES  are 
similar  to  those  for  mean  flow  calculations,  except  at  the  jet  exit.  Here,  the  velocity  field 
must  approximate  a  turbulent  pipe  flow.  Two  attempts  were  made  to  obtain  a  reasonable 
turbulent  inflow. 

In  the  first  effort,  an  auxiliary  calculation  was  performed  of  homogeneous,  isotropic 
turbulence.  This  solution  was  stored  and  used  to  provide  turbulence  which  was  added  to 
the  mean  inflow  velocity  field.  This  approach  should  have  given  a  reasonable  boundary 
condition  for  jet  exit  flow,  but  the  simulation  in  which  it  was  used  was  unsuccessful.  The 
major  shortcomings  in  that  simulation  are  believed  to  be  due  to  factors  other  than  the 
boundary  conditions.  The  simulation  was  run  early  in  the  contract  effort  with  fourth  order 
spatial  discretization  (the  eighth  order  scheme  had  not  yet  been  developed),  and  the  grid 
may  have  been  too  coarse. 

The  second  effort  used  a  random  number  generator  to  specify  the  turbulent 
fluctuations: 


u'  =  Uq  RANF  f(r) 

v'  =  Uq  RANF  f(r)  (23) 

w'  =  Uq  RANF  f(r) 

in  which  u0  was  a  turbulent  scaling  velocity,  typically  about  10%  of  the  jet  maximum 
velocity,  RANF  is  the  Cray  random  number  generator,  and  f(r)  was  a  function  with 
maximum  value  of  unity  which  confines  the  fluctuations  to  the  jet's  shear  layer.  These 
turbulent  fluctuations  were  added  to  the  mean  velocity  and  the  combined  velocity  was 
imposed  at  the  nozzle  exit.  This  form  of  "inflow  turbulence"  is  entire  non-physical,  and  it 
contains  no  inherent  structures.  However,  it  contains  turbulence  energy  in  all  three  velocity 
components,  and  it  provides  a  "seed"  from  which  natural  turbulence  can  evolve. 

Note  that  imposed  perturbations  at  the  nozzle  exit  are  not  required  to  generate  an 
unsteady  simulation.  In  a  study  by  one  of  the  present  authors/3  steady  nozzle  exit 
conditions  were  used  for  a  subsonic  impinging  jet  simulation,  and  that  flow  became 
unsteady.  The  acoustical  energy  emanating  from  the  impingement  region  is  sufficient  to 
perturb  the  jet  shear  layer  at  the  nozzle  exit  and  produce  a  turbulent  flow.  However,  with 
the  limited  resolution  employed  in  three-dimensional  simulation,  the  secondary  instabilities 
which  lead  to  three-dimensional  turbulence  may  be  artificially  suppressed  by  numerical 
errors.  Thus,  the  simulated  turbulence  which  occurs  without  random  forcing  is  prone  to  be 
less  realistic  than  that  with  the  forcing. 

The  validity  of  using  random  forcing  to  replace  the  inflow  turbulence  is  subject  to 
question.  However,  the  results  seem  to  suggest  that  they  are  adequate. 

Time  Integration  -  A  five-stage  second  order  accurate  Runge-Kutta  scheme 
developed  by  Mavriplis  and  Jameson’4  is  used  to  advance  the  equations  in  time.  This 
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scheme  is  explicit  but  is  relatively  stable  and  permits  the  use  of  CFL  numbers  up  to  -2.8 
with  the  present  space  discretization.  The  time  step  can  be  spatially  variable  if  a  calculation 
is  being  run  to  a  steady  state,  as  in  the  case  of  Reynolds-averaged  calculations.  When  time 
accuracy  is  required  for  large  eddy  simulations,  the  time  step  is  fixed,  and  the  maximum 
CFL  number  is  typically  less  than  unity  to  provide  good  temporal  accuracy. 
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4.  LARGE  EDDY  SIMULATIONS 


The  physical  processes  which  control  turbulence  in  impinging  jets  are  poorly 
understood,  as  noted  in  the  introduction.  The  best  experimental  measurement  techniques 
currently  available  cannot  provide  the  information  needed  to  identify  these  processes. 
Hence,  an  alternate  approach  is  required,  and  large  eddy  simulation  (LES)  is  used  in  this 
study.  As  its  name  implies,  LES  involves  flow  simulations  which  resolve  the  large  scales  of 
turbulence  in  addition  to  the  mean  flow.  LES  can  be  used  to  obtain  reasonable 
approximations  of  the  dominant  physical  processes  in  complex  flows  such  as  impinging  jets. 

LES  employs  the  Navier-Stokes  equations  as  the  model  of  the  resolved  portion  of  the 
flow.  The  eighth  order  spatial  scheme  and  the  Runge-Kutta  time  integration,  described 
previously,  are  used  to  integrate  these  equations. 

The  small  scales  of  turbulence  are  not  resolved  in  LES,  and  they  must  be  approximated 
with  some  form  of  model.  Numerous  sub-grid-scale  (SGS)  turbulence  models  have  been 
proposed,  and  they  are  generally  demonstrated  in  flow  such  as  boundary  layers  or 
homogeneous  turbulent  flows.  To  the  present  authors'  knowledge,  no  SGS  model  has  been 
validated  for  turbulence  similar  to  that  in  impinging  jet  flows.  This  is  a  significant  point 
because  the  dominant  turbulent  structure  in  the  impinging  jet  region  consists  of  quasi- 
axisymmetric  vortical  structures.12  It  is  important,  therefore,  that  the  SGS  model  accurately 
predict  the  stresses  in  these  vortices.  Unfortunately,  simple  models  are  notoriously 
inaccurate  for  complex  vortical  flows,  and  typical  SGS  models  will  be  poor  for  LES  of 
impinging  jet  flows. 

No  explicit  SGS  model  is  employed  in  the  present  large  eddy  simulations.  Instead,  the 
artificial  dissipation  operator  produces  diffusive  transport,  which  mimics  small  scale 
turbulent  mixing,  at  the  highest  resolvable  wave  numbers.  The  artificial  dissipation  is  fairly 
specific  to  the  highest  wave  numbers  that  can  be  supported  on  the  grid;  this  permits  the 
simulation  to  capture  the  widest  possible  range  of  turbulent  scales.  The  larger  resolved 
scales  of  turbulence  are  unaffected  by  the  artificial  dissipation. 

LES  Cases  -  Two  flows  were  studied  with  LES;  a  normal  impinging  jet,  and  an 
inclined  impinging  jet  at  conditions  similar  to  the  primary  experimental  test  case.  In  the 
normal  impinging  jet,  the  region  of  "interesting"  physics  is  confined  to  the  impingement 
zone  and  the  initial  part  of  the  radial  wall  jet.  Since  this  region  is  relatively  compact,  it  is 
possible  to  use  a  grid  which  gives  good  numerical  resolution  in  this  region,  while  the  grid 
can  be  coarser  toward  the  outer  boundaries.  In  calculations  of  the  inclined  impinging  jet, 
the  grid  must  have  roughly  uniform  grid  spacing  throughout  the  entire  ground  vortex, 
which  is  a  large  region.  This  requirement  significantly  increases  the  computational  cost  and 
difficulty,  and  reduces  the  accuracy  which  can  realistically  be  achieved.  Hence,  the  case  of 
the  normal  impinging  jet  provides  a  better  understanding  of  turbulence  physics  than  the 
inclined  jet  case.  The  simulation  of  the  inclined  jet  does  reveal  differences  between  the  two 
flows  which  are  critical  for  turbulence  modeling.  The  operating  conditions  of  the  simulated 
inclined  jet  differs  from  the  experiment  The  simulated  Reynolds  number  is  lower  because 
the  grid  resolution  which  can  realistically  be  afforded  is  inadequate  for  a  higher  Reynolds 
number.  The  jet  exit  Mach  number,  Mj,  cannot  be  made  as  small  as  the  experiment  because 
the  stability  limit  on  the  computational  time  step  makes  low  Mach  number  simulations  very 
expensive.  Jet  Mach  numbers  in  the  range  0.4  <  M^  <  0.5  are  low  enough  so  that 
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compressibility  and  aero-acoustic  effects  on  the  turbulence  will  be  negligible.  Neither 
difference  should  materially  affect  the  conclusions  obtain  from  the  simulation  results. 

The  configuration  for  the  normal  jet  LES  is  Hj/Dj  =  3,  Mj  =  0.47  and  ReD  =  50,000.  The 
jet-exit  boundary  layer  thickness  is  6/Dj  =  0.1.  The  inclined  jet  case  is  run  at  H/D=  3.1  and 
l  =  45°,  which  is  the  same  as  the  experiment.  However,  ReD  =  50,000  and  M:  =  0.5.  The  jet 
exit  flow  has  the  same  boundary  layer  thickness  as  the  experiment,  and  tne  turbulence 
energy  and  dissipation  rate  are  appropriate  for  a  turbulent  pipe  flow. 

Dynamical  LES  Results  for  Normal  Impingement  -  LES  generates  large  quantities  of 
data,  and,  at  present,  only  a  fraction  of  that  data  can  be  processed  or  included  in  a  report. 
The  results  presented  here  have  the  specific  purpose  of  illustrating  the  physical  processes 
which  are  critical  to  turbulence  and  which  must  be  modeled  accurately. 

The  dynamical  results  are  presented  first  because  they  provide  good  insight  into  the 
physics  of  turbulence.  Mach  number  contours  on  a  plane  through  the  centerline  of  the  jet 
are  given  in  Figure  12-a.  The  jet  exits  the  nozzle  with  a  relatively  thin  shear  layer  and 
without  any  large  eddy  structure.  At  one-half  the  distance  to  the  impingement  plate,  the  jet 
shear  layer  has  developed  significant  structure,  which  increases  down  to  the  impingement 
plate  and  in  the  wall  jets.  The  maximum  Mach  number  does  not  occur  in  the  free  jet,  but 
rather  in  the  wall  jet  near  x  =  3.5.  This  is  a  somewhat  surprising  result,  but  readily 
explained  by  the  dynamics  of  turbulence 

Contours  of  vorticity  magnitude  on  the  same  plane  are  given  in  Figure  12-b.  The 
relatively  uniform  shear  layer  near  the  nozzle  exit  appears  as  a  continuous  sheet  of  vorticity 
which  then  forms  discrete  vortices.  The  vortices  pair  before  they  strike  the  ground. 
Vortices  in  the  wall  jet  are  seen  at  x  =  ±3.5.  These  are  smaller  and  have  higher  levels  of 
vorticity  than  the  vortices  in  the  jet  above  the  impingement  plate,  despite  the  fact  that  they 
are  older  and  would  therefore  be  larger  and  weaker  due  to  the  diffusion  of  vorticity.  Some 
process  causes  the  vortices  to  intensify  in  the  wall  jet. 

This  impinging  jet  flow  is  axisymmetric  in  the  mean.  Turbulence  in  the  jet  consists  of 
filaments  of  vorticity.  These  filaments  are  not  randomly  oriented,  but  rather  tend  to  form 
rings  that  enclose  the  jet.  Hence,  they  are  quasi-axisymmetric  in  the  free  jet  and  tend  to 
preserve  that  structure  for  some  distance  into  the  wall  jet.  Also,  they  convect  with  roughly 
half  of  the  local  jet  velocity.  When  these  vortex  filaments  flow  outward  in  the  wall  jet,  their 
major  diameter  increases  at  a  rate  roughly  equal  to  Uj/2.  By  conservation  of  mass,  this 
stretching  of  the  major  diameter  causes  the  minor  diameter,  the  thickness  of  the  vortex 
filament  seen  in  the  figure,  to  decrease.  The  vorticity  must  also  increase  to  satisfy  the 
conservation  of  angular  momentum.  A  time  series  of  vorticity  contours  given  in  Figure  13 
indicates  that  the  peak  vorticity  may  triple  as  the  vortices  flow  outward  in  the  wall  jet. 
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Figure  12.  Contours  of  Mach  number  and  vorticity  magnitude  for  normal  impinging  jet. 
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Figure  13.  Time  sequence  of  vorticity  contours  for  normal  impinging  jet. 
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In  addition  to  increasing  the  peak  vorticity,  vortex  stretching  can  increase  the  shear 
stress.  The  instantaneous  velocity  fluctuations  can  be  computed  by  subtracting  the  mean 
velocities  from  the  instantaneous  velocities;  the  instantaneous  turbulent  stresses  are  then 
given  by  the  product  of  two  fluctuating  velocities.  Figure  14  gives  the  instantaneous 
vorticity,  turbulence  energy  (iT iT  +  v'v'  +  w  w  ) / 2.  and  the  shear  stress,  u  v  .  The  peaks  of 
energy  and  shear  stress  are  associated  with  the  vortices:  as  expected,  the  turbulent  stresses 
are  carried  by  the  vortices.  The  maximum  values  of  energy  and  shear  stress  at  this  instant 
are  larger  in  the  wall  jet  than  in  the  free  jet.  This  trend  also  applies  to  the  time-averaged 
shear  stress  as  shown  in  statistical  results  presented  below. 
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(b)  Turbulence  energy 
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(c)  Shear  stress 

Figure  14.  Contour  plots  on  instantaneous  vorticity,  turbulence  energy  (u'u'  + 
v'v'  +  w'w')/2  and  the  shear  stress,  u'v'. 


Hence,  vortex  stretching  is  a  critical  mechanism  which  affects  turbulence  in  the 
impingement  region,  particularly  in  the  wall  jet  flowing  away  from  the  impingement  point. 

It  is  reasonable  to  question  whether  the  above  implication  of  vortex  stretching  as  a 
critical  turbulence  mechanism  results  from  approximation  used  in  the  numerical  simulation 
method.  The  importance  of  vortex  stretching  would  be  over-emphasized  in  a  simulation 
which  is  axisymmetric.  Therefore,  it  is  important  to  establish  the  three-dimensionality  of 
the  flow,  both  in  the  free  jet  and  in  the  wall  jet.  Figures  12  and  13  show  that  the  free  jet  is 
not  symmetric  about  its  centerline.  Figure  15  give  the  instantaneous  Mach  number  contours 
on  a  plane  normal  to  the  free  jet,  and  again,  significant  highly  three-dimensional  structures 
can  be  seen.  A  top  view  of  Mach  number  contours  and  velocity  vectors  in  the  wall  jet 
(y  =  0.3)  in  Figure  16  also  shows  that  the  flow  is  very  much  three-dimensional;  however,  a 
dominant  quasi-axisymmetric  structure  does  exist.  The  velocity  vectors  suggest  that  the 
three-dimensionality  increases  with  increasing  radius.  These  results  demonstrate  that  three- 
dimensional  turbulence  exists  in  the  simulations.  This  suggests,  but  does  not  prove,  that  the 
vortex  stretching  mechanism  is  not  artificially  augmented  by  the  limited  resolution  available 
in  the  simulation. 

Dynamical  LES  Results  for  the  Primary  Test  Case  -  A  large  eddy  simulation  was  also 
performed  of  the  primary  test  case.  As  mentioned  above,  the  magnitude  of  this  calculation 
was  significantly  greater  than  the  normal  impinging  jet.  The  calculation  was  performed  on 
a  grid  which  had  less  resolution  than  the  previous  case,  and  the  time  scale  of  the  global  flow 
field  was  significantly  longer  than  the  flow  in  just  the  impingement  region.  This  calculation 
was  not  run  long  enough  to  obtain  statistics.  However,  information  concerning  the 
dynamical  structure  of  turbulence  was  obtained,  which  provides  insight  into  the  important 
vortex  stretching  process. 
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Figure  15.  Instantaneous  Mach  number  contours  on  a  plane  normal  to  the  jet  and 
above  the  impingement  zone,  for  the  normal  impinging  jet. 
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(b)  Velocity  vectors 


Top  view  of  instantaneous  Mach  number  contours  and  velocity  vectors 
at  y  =  0.3,  for  normal  impinging  jet  case. 


(b)  Vorticity  magnitude,  time  =•  55 
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(c)  Vorticity  magnitude  =  64 

Figure  17.  Contours  of  Mach  number  and  vorticity  magnitude  on  the  center  plane 
for  inclined  jet  case,  and  voriticity  contours  at  a  later  time  showing 
vortex  evolution. 

Figure  17  gives  Mach  number  and  vorticity  magnitude  contours  on  the  centerplane 
through  the  jet  pipe  and  ground  vortex  at  one  instant  in  time.  The  Mach  number  contours 
reveal  the  general  location  of  the  jet  shear  layers,  and  the  fact  that  the  majority  of  the  jet 
fluid  goes  forward  into  the  ground  vortex.  The  dynamical  structure  due  to  turbulence  is  not 
as  pronounced  as  in  the  normal  impinging  jet  case,  perhaps  because  of  the  reduced  grid 
resolution.  As  in  the  case  of  the  normal  impinging  jet,  vortex  stretching  amplifies  the  peak 
vorticity  and  decreases  the  vortex  size  in  the  wall  jet.  Vorticity  contours  in  Figure  17-b  and, 
at  a  later  time  in  Figure  17-c,  reveal  several  discrete  vortices.  The  vortices  upstream  of  the 
impingement  point  (x  <  -5.0)  are  comparable  in  strength  to,  or  weaker  than,  the  peak 
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vorticity  in  the  jet  shear  layer.  Comparison  of  Figures  17-b  and  17-c  reveals  that  vortices  do 
intensify  in  this  region;  however,  the  peak  vorticity  is  significantly  less  amplified  than  in  the 
normal  impingement  case.  Downstream  of  the  impingement  point  several  compact,  intense 
vortices  can  be  seen.  The  peak  turbulent  vorticity  in  Figure  17-c  occurs  in  the  wall  jet  at  x  =  - 
2.5.  (This  figure  also  reveals  the  scale  of  vortices  that  can  be  resolved  on  this  grid.)  A  top 
view  of  the  velocity  vector  field  is  given  in  Figure  18.  The  lateral  divergence,  -dw/dz, 
appears  to  be  stronger  downstream  of  the  impingement  point  than  upstream  of  it.  This 
characteristic  of  the  mean  flow  causes  the  different  amplification  of  vortices  upstream  and 
downstream  of  the  impingement  zone. 
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Figure  18.  Velocity  vectors  at  y  =  0.03  (at  every  other  grid  point)  showing  lateral 
divergence,  -3w/dz,  in  regions  of  vortex  intensification. 


Statistical  Results  for  Normal  Impingement  -  In  order  to  obtain  turbulent  statistics, 
the  turbulence  simulation  was  run  until  it  appeared  to  have  reached  an  asymptotic  unsteady 
condition,  independent  of  its  initial  conditions.  Integration  was  then  started  to  obtain  time- 
averages  of  various  flow  variables.  This  integration  was  performed  over  a  reasonable 
length  of  time,  until  the  turbulent  statistics  computed  from  these  averages  demonstrated 
clear  trends.  In  these  simulations,  the  integration  was  not  carried  out  long  enough  to  obtain 
converged  statistics  because  the  computational  cost  would  have  exceeded  the  available  CPU 
time.  For  this  flow  which  is  steady  and  axisymmetric  in  the  mean,  the  turbulent  statistics 
should  display  an  appropriate  symmetry  about  the  jet  centerline.  As  will  be  seen  below,  the 
results  are  not  symmetrical  because  the  statistics  are  not  converged.  However,  the 
integration  time  was  sufficient  to  obtain  an  understanding  of  the  dominant  physical 
characteristics  of  impinging  jet  turbulence,  and  to  identify  some  critical  differences  between 
normal  and  inclined  impinging  jets. 

These  impinging  jet  flows  are  three-dimensional  in  the  mean.  Averages  are  computed 
with  respect  to  time;  there  are  no  spatial  directions  in  which  the  mean  flow  or  turbulence  is 
homogeneous  so  that  spatial  averaging  can  be  used  in  addition  to  temporal  averaging.  The 
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calculation  of  turbulent  statistics  requires  a  large  amount  of  memory  and  integration  time. 
To  reduce  the  memory  requirements,  the  statistics  are  computed  at  every  other  grid  point  in 
each  spatial  direction;  this  reduces  the  memory  required  for  the  statistics  calculation  by  a 
factor  of  8.  Time  averages  were  obtained  from  the  simulation  which  permitted  the 
calculation  of  the  following  statistics: 

p,  Uj ,  p  mean  variables 

utUj  turbulent  stresses 

pduj  /  dXj  pressure  strain  correlation 

which  are  believed  to  be  the  quantities  of  greatest  importance  in  these  flows. 

The  simplest  statistics  are  the  mean  flow  quantities,  and  the  mean  Mach  number  and 
pressure  contours  and  the  velocity  vectors  are  given  in  Figure  19.  Qualitatively,  these 
results  are  similar  to  results  from  a  conventional  Reynolds-averaged  calculation 


(a)  Mach  number 
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(b)  Pressure  contours 


(c)  Velocity  vectors 

Figure  19.  Mean  Mach  number  and  pressure  contours,  and  velocity  vectors  on  center 
plane  of  normal  impinging  jet,  computed  from  time  average  of  LES  results. 


Contours  of  turbulence  kinetic  energy  on  a  plane  through  the  center  of  the  jet  are  given 
in  Figure  20-a.  The  turbulence  energy  is  below  the  minimum  contour  level  at  the  nozzle  exit 
but  increases  rapidly  in  the  shear  layers  of  the  free  jet.  The  peak  energy  is  reached  in  the 
impingement  zone  and  remains  relatively  large  in  the  wall  jet.  The  u'u'  and  v'v'  normal 
stresses  in  Figure  20-b  and  20-c  reveal  significantly  greater  spatial  variations  than  the 
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turbulence  energy  in  the  impingement  zone.  The  u'u'  (which  is  parallel  to  the  wall)  has 
local  maxima  very  near  the  wall,  beneath  the  wall  jet,  and  in  the  impingement  zone.  The 
maximum  v'  v'  values  also  occur  in  the  wall  jet,  but  somewhat  above  the  wall  at  a  location 
that  corresponds  roughly  with  the  vortex  cores  seen  in  Figures  12  and  13.  The  shear  stress, 
given  in  Figure  20-d,  also  reaches  its  maximum  values  in  the  wall  jet  at  about  the  same 
locations  as  the  maximum  in  v'v'.  These  results  further  emphasize  the  relation  between 
unsteady  vorticity  and  the  turbulent  stresses,  and  hence,  the  significance  of  vortex 
stretching  as  a  means  of  increasing  shear  stress  in  the  wall  jet. 


(b)  Normal  stress,  u'  u 


oiiii'iiaiii 

5  —0.0  —4.0  —3.0  -2.0  —  1.0.00*10**0.0  2.0  3.0  4.0  0.0 

o 

X 


(c)  Normal  stress,  v'  v' 
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(d)  Shear  stress,  u'  v' 


Figure  20.  Contour  plots  of  turbulent  statistics  on  center  plane  of  normal  impringing  jet, 
computed  from  time  average  of  LES  results. 


There  is  also  a  tendency  for  the  shear  stress  to  increase  in  the  wall  jet  while  the 
turbulence  energy  does  not.  The  correlation  coefficient  R12  =  u'v'  /(u'u'  v'v')Vj  is  also 
larger  in  the  wall  jets  than  in  the  free  jet.  These  results  indicate  that  some  process  in  the  wall 
jet,  presumably  vortex  stretching,  causes  turbulence  there  to  be  more  efficient  at  producing 
shear  stress  from  the  available  turbulence  energy.  This  has  important  implications  for 
turbulence  modeling,  because  it  implies  that  vortex  stretching  can  cause  the  shear  stress  to 
increase  without  increasing  the  turbulence  energy.  Vortex  stretching  modifications  to  the 
k-e  model  should  therefore  be  directed  principally  at  the  e— equation.  However, 
consideration  of  the  Reynolds-stress  transport  equations  indicate  that  this  is  not  strictly  so. 

In  summary  of  this  section,  the  LES  results  have  provided  dynamical  and  statistical 
information  about  impinging  jet  turbulence.  Both  viewpoints  indicate  that  vortex  stretching 
is  a  critical  process  in  the  wall  jet  which  emanates  from  the  impingement  region.  Vortex 
stretching  intensifies  vorticity  fluctuations,  and  certain  components  of  turbulence  are  also 
increased.  A  significant  difference  between  the  normal  and  inclined  impinging  jets  is  the 
magnitude  of  vortex  stretching  in  the  wall  jet  that  flows  forward  (into  the  crossflow)  from 
the  impingement  point.  The  wall  jet  from  the  inclined  jet  experiences  less  lateral  divergence, 
and  the  shear  stress  is  less  amplified,  than  in  the  normal  impinging  jet. 
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5.  TURBULENCE  MODELING 


The  k-e  model  is  used  for  the  basis  of  this  work  because  it  is  one  of  the  best  all-purpose 
models  currently  available,  because  it  is  relatively  easy  to  use,  and  because  more  complex 
models  have  not  been  demonstrated  to  be  superior  for  impinging  jet  flows.  As  noted  in  the 
Introduction,  many  types  of  turbulence  models  have  been  applied  to  impinging  jet  flows, 
and,  almost  without  exception,  these  models  have  failed  to  give  adequate  accuracy. 
Modifications  to  the  k— e  model  which  account  for  the  effects  of  streamline  curvature  and 
vortex  stretching  are  described  in  this  section.  These  modification  should  also  be  applicable 
to  many  other  two-equation  eddy  viscosity  models  via  a  transformation  of  variables. 

The  Standard  k— e  Model 

Turbulent  stresses  are  computed  via  the  Boussinesq  approximation. 


pujuj  =-jit<dU1/dx,+  dU,/dx1)+  6ij2pk/3  (24) 


in  which 


Ht  =  Cllk2/e 


(25) 


is  the  eddy  viscosity.  The  k-equation  is  obtained  from  the  trace  of  the  Reynolds-stress 
transport  equations  with  an  approximation  for  diffusive  transport  and  is  believed  to  be 
without  major  flaws.  It  can  be  written  as 


Dpk 

Dt 


ixj  "k  ax3 


(26) 


in  which  Pk  is  the  production  of  turbulence  energy  derived  by  contracting  the  production 
tensor  in  the  Reynolds-stress  transport  equations 


pij  =  -  P  “iHn  *Uj/dxn,  -  P  UX 


The  e~ equation  is 


Pk  =  -p  Uj Uj  (dU/dXj  +  dUj/aXj)/2 


_  „  d  u  ds 

Dpg _ _t  _  _  g 

DC  dx.  a  dx.  E 

3  g  j 


in  which  the  source  term  for  the  standard  model  is  given  by 


S  =  C  p—  P  -  C 
e  erk  k  z2  k 


(27-a) 

(27-b) 


(28) 


(29) 
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The  e-equation  was  originally  derived  from  the  fluctuating  vorticity  transport  equation.  It 
can  alsr  be  derived  from  the  Reynolds-stress  transport  equations  for  shear  stress  and 
energy,  using  the  k-e  expression  for  shear  stress,  Eq.  (24),  to  obtain  the  t  transport 
equation.  Equations  (24-29)  define  the  standard,  high-Reynolds  number  form  of  the 
k—  e  model.  Model  coefficients  are  given  in  Table  1.  These  equations,  with  modifications 
described  below,  are  solved  throughout  the  flow  domain  except  at  the  wall. 

Between  the  point  at  the  wall  and  the  first  interior  point,  the  so-called  wall  functions 
are  used  to  approximate  the  mean  flow  and  turbulence.  With  the  assumptions  of 
incompressible  flow,  a  logarithmic  velocity  profile,  constant  shear  stress,  and  turbulence 
production  equaling  dissipation,  the  wall  function  relations  can  be  derived.  (See,  e  g..  Ref. 
15.)  These  equations  give  wall  boundary  conditions  for  k  and  e,  and  the  wall  shear  stress. 

The  standard  k—  e  model  works  reasonably  well  for  a  limited  range  of  flows,  but  its 
accuracy  is  poor  for  most  complex  flows.  Modifications  which  improve  the  model's 
accuracy  are  described  below. 

Modifications  to  the  k—  i  Model 

An  underlying  tenet  of  the  present  work  is  that  the  effects  of  streamline  curvature  and 
vortex  stretching  are  two  of  the  most  significant  extra  strains  which  occur  in  impinging  jet 
flows,  and  which  are  not  treated  accurately  by  the  standard  k-e  model.  If  the  turbulence 
model  can  be  modified  so  that  it  performs  well  for  these  phenomena,  for  fundamentally 
correct  reasons,  it  should  also  perform  well  in  complex  impinging  jet  flows. 

Streamline  Curvature  -  There  is  a  considerable  body  of  research  dealing  with  the 
effects  of  streamline  curvature  on  turbulence.16-20  The  present  model  is  similar  to  the  one 
proposed  by  Launder  et  al.18  At  the  heart  of  the  curvature  model  is  the  curvature  stability 
parameter  p  which  indicates  whether  the  additional  strains  due  to  streamline  curvature 
promote  or  suppress  turbulence.  For  a  simple  curved  flow 

P  =  -(2V/R)(aV/dn  +  V/R)  (30) 

in  which  V  is  the  mean  fluid  speed,  R  is  the  streamline  radius  of  curvature,  and  n  is  the 
outward  normal  to  the  streamline. 

It  is  beneficial  to  recast  the  stability  parameter  in  a  form  which  is  easily  computed  in  a 
general  three-dimensional  flow.  The  expression  for  p  can  be  written  as 

0  =  -(2/R)(V2/R  +  v2aV2/dn)  (31) 

The  pressure  gradient  normal  to  streamlines  in  a  steady  curved  flow  is  given  by 

8p/an  =  pV2/R  (32) 

so  that  Eq.  (31 )  becomes 

P  =  -(2/pr)  dp,/ dn  (33) 
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in  which  p,  =  p  +  pV2/ 2  is  the  incompressible  total  pressure.  Equation  (33)  reveals  that 
potential  flow,  for  which  pt  is  constant,  is  neutrally  stable.  There  are  at  least  two  possible 
ways  to  compute  0.  One  is  to  determine  the  streamline  curvature  from  the  pressure 
gradient,  as  in  Eq.  (32),  which  yields 

P  =  -2  dp/ 8n  dpt/dn  /(jfiV2)  (34) 

In  flows  where  viscous  effects  are  weak,  the  total  pressure  varies  slowly  along  streamlines 
so  that  Eq.  (34)  can  be  approximated  by 

P  =  -2  dp/dXj  dpf/dXj  /( p 2V2)  (35) 

Streamline  curvature  1  /  R  may  also  be  computed  from  the  flow  geometry  as 

l/R=dG/ds  (36) 

in  which  6  is  the  flow  angle  and  s  the  distance  along  a  streamline.  In  three  dimensions 
curvature  has  both  magnitude  and  orientation  and  is  thus  a  vector.  The  curvature  vector 
pointing  from  the  center  of  curvature  outward  can  be  expressed  as 

i/Ri  =  -uj/va(Ui/v)/axj  (37) 


in  which  V  is  the  velocity  magnitude,  so  that  Uj/V  is  the  flow  direction.  Taking  the  inner 
product  with  dpt/dx{  as  in  Eq.  (35)  yields 

p  =  (2/p)  (U/V)  aoyvvaxj  dp^  (38) 

Either  Eq.  (35)  or  (38)  can  be  used  to  compute  p,  and  they  give  similar  values  except 
where  the  assumptions  used  to  derive  them  are  invalid,  such  as  through  shocks,  in  other 
highly  viscous  regions,  or  where  the  velocity  goes  to  zero.  Equation  (38)  tends  to  be  better 
behaved  in  general  because  it  minimizes  the  problem  of  the  denominator  going  to  zero  at 
stagnation  points,  and  it  is  used  in  the  present  calculations. 

The  stability  parameter  is  then  used  to  modify  the  £ -equation,  by  adding  the  term 

Cckp  (39) 

to  the  source  term  S£.  This  choice  for  the  form  of  the  additional  source  term  is  not  arbitrary, 
but  has  been  found  to  perform  better  than  another  dimensionally  correct  form, 
Cc  t  I p I v*  sign {p). 

Vortex  Stretching  (Lateral  Divergence)  -  Vortex  stretching  is  often  referred  to  as 
lateral  divergence  because  it  typically  involves  the  lateral  stretching  (dW/dz  >  0)  of  a  two- 
dimensional  shear  layer  in  the  x-y  plane.  By  conservation  of  mass  and  angular  momentum, 
vortex  stretching  intensifies  vorticity.  The  effects  of  vortex  stretching  on  the  mean  flow  are 
easily  treated  with  the  Euler  equation.  However,  its  effects  on  turbulence  are  complex  and 
not  well  understood.  Obviously,  therefore,  modeling  for  the  vortex  stretching  process  is  not 
well-developed. 
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Some  first  order  effects  of  vortex  stretching  can  be  deduced  from  Eq.  {27),  the 
Reynolds-stress  production  term.  For  the  above  example  of  a  typical  shear  layer.  Pi. 
indicates  that  vortex  stretching  causes  a  transfer  of  turbulence  energy  from  the  w'  w  7 
component  into  the  u'u'  and  v' v'  components,  and  that  the  details  of  this  transfer  depend 
upon  the  initial  turbulence  and  the  strain  field.  In  this  shear  flow,  the  increase  inv'v' 
interacts  with  the  mean  shear  to  increase  the  production  of  the  shear  stress.  The  anisotropic 
normal  stresses  interacting  with  the  mean  strain  cause  vortex  stretching  to  transfer  energy 
from  the  mean  flow  to  the  turbulence,  especially  for  strong  vortex  stretching.  Hence,  vortex 
stretching  causes  the  turbulence  energy  and  shear  stress  to  be  higher  than  in  a  similar  flow, 
but  with  8W/dz  =  0.  Consideration  of  a  single  vortex  subject  to  stretching  leads  to  similar 
conclusions. 

The  effect  of  vortex  stretching  should  depend  significantly  on  the  state  of  the 
turbulence,  for  example,  the  distribution  of  the  energy  among  the  normal  stresses  and  the 
structure  of  the  vortices.  Wall-bounded  flows,  round  and  plane  jets,  and  low-  and  high¬ 
speed  free  shear  layers  have  very  different  vortical  structures  and,  hence,  should  respond 
differently  to  vortex  stretching. 

Modeling  for  the  effects  of  vortex  stretching  is  in  a  primitive  state,  perhaps  because  of 
the  complexities  mentioned  above.  Pope21  suggested  vortex  stretching  as  the  cause  of  the 
"round  jet/plane  jet"  anomaly,  in  which  the  round  jet  spreads  less  rapidly  than  the  plane  jet, 
but  the  standard  k-e  model  predicts  the  opposite  trend.  The  proposed  vortex  stretching 
model  was  an  additional  source  term  in  the  e-equation 

C^GjtkVe  (40) 

in  which  the  vortex  stretching  is  given  by, 

x  =  n,njsij  (41) 

the  vorticity  is 

n^e^aiyax^  (42) 


with  e^  being  the  alternating  third  rank  tensor,  and  the  strain  rate  tensor  is  given  by 

Stj  =  (dUi/dXj+dU/axjVZ  (43) 

In  Pope's  model  the  function  G  =  1.  The  model  coefficient  was  optimized  to  be  Cvs  =  0.79, 
and  the  model  then  gave  the  correct  spreading  rate  for  the  round  and  plane  jets.  In  a  flow  in 
which  vortex  stretching  occurs,  Eq.  (40)  increases  the  dissipation  which  reduces  the  eddy 
viscosity. 

Unfortunately,  this  model  produces  the  wrong  trend  in  a  radial  wall  jet,22  which  is  a 
flow  of  major  interest  for  impinging  jet  predictions,  and  in  the  radial  free  jet.  The  shear 
stress  in  these  two  flows  is  considerably  higher  than  in  their  two-dimensional  counterparts, 
whereas  the  Pope  modification  reduces  the  modeled  shear  stress. 
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To  make  the  k—  e  model  sensitive  to  the  vortex  stretching  which  is  so  pervasive  in 
impinging  jet  flows,  a  modification  of  the  form  of  Eq.  (40)  must  be  used.  Note  that  the  sign 
of  GyS  required  for  impinging  jets  is  opposite  of  that  used  by  Pope.  A  further  difficulty  must 
also  be  corrected.  Apparently  because  of  their  differing  vortical  structure,  the  free  and  wall- 
bounded  jets  require  slightly  different  modeling.  An  ad  hoc  remedy  is  to  make  G  a  function 
of  wall  proximity. 


G(L/y)  =  1.  +  Cvs2L  /  y  (44) 

in  which  L  =  C<l3/4k3/2/c  is  the  turbulent  mixing  length  and  y  is  the  distance  from  the  wall. 
To  ensure  the  model  is  robust,  L  is  limited  to  L  s  y.  This  form  of  the  model  makes  it 
possible  to  match  the  growth  rates  of  free  and  wall-bounded  radial  jets.  A  less  empirical 
wall  treatment  would  require  knowledge  of  the  turbulence  structure  which  is  not  readily 
available  from  the  k—  e  model. 

The  full  source  term  for  the  e-equation,  including  curvature  and  vortex  stretching 
effects,  can  be  written  as 


=  Cel^k 


Pk  - 


,  ££_ 

'e2  k 


C  k  p 
c  r 


♦  CvsG(L/y>,  ~ 


(45) 


Disadvantages  of  Vortex  Stretching  Model  -  Despite  the  need  for  a  model  of  vortex 
stretching  effects  on  turbulence,  the  model  presented  above  can  cause  difficulties  in 
calculations  of  some  flows.  As  noted  previously,  the  physics  of  vortex  stretching  is  complex 
and  poorly  understood,  and  only  a  few  efforts  have  specifically  addressed  the  issue  of 
modeling  for  this  phenomenon. 

For  radial  jets,  Cvs  must  have  the  opposite  sign  of  that  used  by  Pope.21  However,  this 
value  of  Cvs  exacerbates  the  problem  which  Pope  addressed:  the  k—  e  model's  error  in  a 
round  jet.  The  wall-proximity  function  F(L/y)  causes  the  vortex  stretching  effects  to  be 
relatively  small  away  from  the  wall,  so  that  this  extra  error  in  the  round  jet  remote  from 
walls  is  also  small.  It  is  now  well-known  that  coherent  structures  dominate  the  momentum 
transport  by  turbulence,  and  that  these  structures  differ  from  flow  to  flow.  For  example,  the 
dominant  structures  in  low-speed  plane  shear  layers  consist  of  spanwise  vortices,  while 
fully  developed  round  jets  contain  helical  vortices.  It  has  been  proposed23  that  helical  or 
swept  vortices  are  less  efficient  than  spanwise  ones  at  generating  mixing.  Thus,  vortex 
orientation  may  be  the  cause  of  the  round  jet-plane  jet  anomaly,  rather  than  vortex 
stretching.  The  k-c  model  has  no  "knowledge"  of  vortical  structure  and  cannot  model  the 
effects  of  structural  changes.  There  remains  a  fundamental  and  practical  issue  concerning 
the  k— £  model's  accuracy  for  round  jets. 

A  second  difficulty  is  related  to  the  mathematically  ill-posed  nature  of  the  k—  £  model. 
Ratios  of  k  and  £  which  appear  many  places  in  this  model,  and,  in  laminar  flow,  these  ratios 
are  singular.  Physically  these  ratios  must  be  reasonably  well-behaved.  However,  numerical 
solutions  of  the  turbulence  model  equations  do  not  provide  any  such  guarantees,  and  many 
codes  which  use  the  k-c  model  employ  artificial  limits  which  insure  that  k,  £  >  0.  These 
limits  do  not  impart  physically  correct  behavior  to  the  model;  they  merely  prevent  the 
modeled  solution  from  diverging,  under  most  circumstances.  The  presence  of  £  in  the 
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denominator  of  Eq.  (40)  gives  the  vortex  stretching  term  the  ability  to  exhibit  ill-conditioned 
behavior.  Vortex  stretching  causes  t  to  decrease  to  provide  a  higher  shear  stress,  but  this 
change  also  increases  the  magnitude  of  the  vortex  stretching  term.  This  nonlinear  decrease 
of  t,  and,  hence,  increase  of  fit,  in  the  presence  of  vortex  stretching  may  be  somewhat 
physical;  analysis  of  the  production  terms  in  the  Reynolds-stress  transport  equations 
indicates  that  vortex  stretching  is  a  nonlinear  process.  However,  the  limits  to  stress 
amplification  which  govern  physical  turbulence  are  not  present  in  the  k-c  model. 

One  of  the  conclusions  of  this  work  is  that,  while  vortex  stretching  is  a  critical  process 
which  must  be  included  in  a  model  of  turbulence,  a  deeper  understanding  of  the  physical 
processes  and  a  mathematically  well-posed  model  should  be  pursued  in  future  basic 
research  in  this  area. 

Model  Nomenclature  -  For  brevity  and  clarity,  the  following  abbreviations  will  be 
used  in  discussing  the  results  given  with  possible  forms  of  the  k—  e  model. 

STD  standard  model  defined  by  Eqs.  (24-29), 

C  with  curvature  modification,  Eqs.  (38, 39), 

VS  with  vortex  stretching  modification,  Eqs.  (40-44), 

CVS  with  curvature  and  vortex  stretching  modifications. 

For  illustrative  purposes,  the  different  possible  combinations  of  model  modifications  will  be 
used  to  compute  flows.  However,  the  CVS  model  is  believed  to  be  the  most  accurate  of 
those  considered  here,  and  it  should  be  used  for  all  impinging  jets  flow  calculations. 

Model  Calibration  -  The  new  model  coefficients  introduced  with  the  curvature  and 
vortex  stretching  modifications  are  determined  by  numerical  optimization  for  simple  flows 
that  are  specific  to  the  phenomena  of  streamline  curvature  and  vortex  stretching. 

The  three  flows  used  in  the  model  calibration  are  self-similar  in  the  streamwise 
direction:  a  wall  jet  over  a  convex  logarithmic  spiral,  a  radial  wall  jet,  and  a  radial  free  jet. 
Consensus  values  for  the  growth  rate  and  maximum  shear  stress  in  wall  jets  were  obtained 
from  a  review  by  Launder  and  Rodi,24  while  the  free  jet  data  is  from  Witze  and  Dwyer.25  A 
computer  program  for  self-similar  flows  was  used  so  that  very  good  numerical  resolution 
was  possible.  The  computed  results  were  grid-independent.  The  standard  model 
coefficients  C^,  Clc,  C2e,  <7k,  ere,  and  Pr,  are  unchanged  from  the  values  used  by  many 
modelers.15  Prt  is  the  turbulent  Prandtl  number.  The  model  coefficients  are  given  in  Table  2. 

Table  2.  Model  Coefficients. 

Cp  Cel  Ce2  Cc  Cvs  Cvs2  <rk  at  Prf 

0.09  1.44  1.92  0.2  -0.15  4.0  1.0  1.3  0.9 

Cc  is  similar  to  the  values  selected  by  Launder  et  al18  and  Rodi  and  Scheuerer20  for  several 
different  flows,  while  Cvs  is  very  different  from  previously  used  values.21,26 

A  comparison  of  experimental  data  and  computed  model  results  for  the  growth  rate  of 
the  curved  wall  jet  is  given  in  Figure  21.  At  zero  curvature  the  growth  rate  dy1;j  /  dx  =  0.073, 
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in  which  yv  is  the  half-velocity  thickness.  The  k-£  model  gives  a  growth  rate  of 
dyiJ} !  dx  =  0.69,  which  is  23%  high  and  which  is  typical  for  this  model.  Cc  was  selected  to 
match  the  change  of  growth  rate  with  respect  to  curvature,  d(8y  XJ  dx)/d(ylf  /R),  and, 
hence,  it  is  not  surprising  that  the  CVS  model  performs  well  over  the  full  range  of 
streamline  curvature.  Also  included  for  comparison  are  results  from  the  STD  model,  which 
gives  only  about  25%  of  the  increase  in  growth  rate  due  to  curvature.  The  slight  increase  in 
the  growth  rate  is  due  to  the  additional  shear  associated  with  the  streamline  curvature. 


Figure  21 .  Growth  rates  of  self-preserving  wall  jet  on  a  logarithmic  spiral,  experiment  and 
calculations  with  the  STD  and  C  models. 


The  measured  and  predicted  growth  rates  of  the  free  and  wall  jets  in  planar  and  radial 
configurations  are  given  in  Table  3.  The  standard  model  underpredicts  the  growth  rate  of 
the  radial  jets  by  a  substantial  amount.  Cvs  and  Cvs2  were  selected  to  reproduce  the  change 
in  growth  rate  due  to  streamline  curvature,  and,  because  of  the  model's  error  for  the  plane 
wall  jet,  the  CVS  model  does  exhibit  some  error  for  these  radial  flows. 

Table  3.  Calibration  Flows  for  Vortex  Stretching  Model 


Growth  Rate,  dyVj/dx 


Exp. 

STD  Model 

CVS  Model 

Plane: 

Free  jet 

0.10-0.11 

0.11 

— 

Wall  Jet 

0.073 

0.09 

- 

Radial: 

Free  Jet 

0.10 

0.087 

0.11 

Wall  Jet 

0.09 

0.54 

0.10 
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6.  REYNOLDS-AVERAGED  IMPINGING  JET  CALCULATIONS 


Results  are  presented  for  Reynolds-averaged  calculations  of  three  impinging  jets  in 
crossflow.  The  experiments  of  these  flows  were  conducted  in  small  wind  tunnels.  In  the 
calculations,  the  tunnel  geometry  has  been  represented  as  faithfully  as  possible.  The  tunnel 
side  walls  are  included  in  all  calculations,  and  the  wind  tunnel  top  wall  is  included  in 
Cases  I  and  II.  For  Case  III,  the  tunnel  section  was  rather  tall,  roughly  ten  times  the  height 
of  the  ground  vortex.  A  shorter  computational  domain  was  used  with  free  inflow-outflow 
boundary  conditions  imposed  on  the  top  of  the  domain  to  reduce  the  computational  cost. 
The  pipe  which  supplies  the  jet  has  been  approximated  with  interior  boundary  conditions  in 
all  cases.  While  the  experiments  are  at  incompressible  speeds,  the  calculations  were  run  at 
intermediate  jet  Mach  numbers,  typically  Mj  =  0.5,  to  maintain  computational  efficiency 
with  the  compressible  algorithm;  compressibility  effects  in  the  computed  results  are 
negligible. 

Case  I:  Normal  Jet,  HID  =  3,  V./V.  =  5 

This  case  is  of  a  jet  with  a  nearly  uniform  velocity  issuing  from  a  nozzle  normal  to  the 
plate.  The  crossflow  causes  significant  deflection  before  the  jet  strikes  the  ground.  This  flow 
has  been  studied  experimentally  by  Cimbala  et  al,12  using  laser  velocimetry,  surface 
pressures  and  flow  visualization. 

The  velocity  magnitude  contours  and  velocity  vectors  from  the  solution  computed  with 
the  CVS  model  are  given  in  Figure  22.  The  jet  flows  downward  from  x/D)  =  0  and  y/Dj  =  3. 
The  cross  flow  from  the  left  in  the  figure  causes  the  ground  vortex  to  form  on  the  left  of  the 
impingement  point.  The  general  features  of  the  flow  are  easily  identified  in  this  figure. 
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(a)  Velocity  magnitude 
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(b)  Velocity  vectors 

Figure  22.  Velocity  magnitude  contours  and  velicity  vectors  on  the  center  plane, 
computed  with  the  CVS  model,  for  Case  I. 

Surface  centerline  pressure  coefficients  from  experiment  and  calculations  are  given  in 
Figure  23.  The  high  pressure  corresponds  to  the  impingement  point,  while  the  low  pressure 
upstream  of  the  impingement  point  indicates  the  position  of  the  ground  vortex.  For  this 
flow,  the  STD  model  predicts  the  ground  vortex  to  be  too  far  upstream  of  the  experimental 
location  by  about  50%.  (Upstream  is  with  respect  to  the  crossflow.)  This  type  of  error, 
including  direction  and  magnitude,  has  also  been  obtained  by  other  researchers,  as  noted  in 
the  Introduction.  The  C  model  predicts  the  ground  vortex  to  be  even  further  upstream  of 
the  experimental  location,  which  increases  the  error  relative  to  the  STD  model.  The  VS 
model  predicts  the  ground  vortex  to  be  in  virtually  the  correct  location,  while  the  CVS 
model  places  the  ground  vortex  upstream  of  the  experimental  location  by  about  13%. 

A  comparison  of  measured  velocity  vectors  and  those  computed  with  the  CVS  model  in 
the  ground  vortex  is  given  in  Figure  24.  As  in  Figure  23,  the  computed  ground  vortex 
location  is  slightly  upstream  of  the  experimental  location,  which  is  consistent  with  the 
surface  static  pressures. 
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(b)  Measured  velocities 


Figure  24.  Comparison  of  center  plant  mean  velocity  vectors  for  Case  I  computed  with 
the  CVS  Model  and  measured  by  Cimbala  et  al.29 


The  surface  pressure  measurements  with  the  VS  model  are  very  close  to  the 
experimental  data,  but  this  model  modification  alone  cannot  be  truly  correct.  Despite  the 
small  degradation  in  accuracy,  the  curvature  modification  must  be  included  in  calculations 
of  this  flow.  The  curvature  modification  is  required  to  compute  accurately  the  shear  stresses 
in  the  impingement  zone  of  a  normally  impinging  jet,  the  curved-wall  jet,  and  other 
curved  flows.  This  result  also  indicates  that  the  apparent  marginal  accuracy  of  the  STD 
model  is,  in  fact,  due  to  luck,  because  the  errors  incurred  by  not  accurately  predicting 
curvature  and  vortex  stretching  effects  are  of  opposite  sign. 

The  significance  of  the  CVS  model  modifications  to  the  computed  turbulence  levels  can 
be  seen  in  the  predicted  centerplane  shear  stress,  given  in  Figure  25.  The  STD  model 
predicts  the  shear  stress  in  the  vortex  core  (x/Dj  =  -2.5,  y/Dj  =  0.5)  to  be  very  nearly  the 
maximum  in  the  domain,  exceeding  the  stress  in  the  free  jet  by  more  than  a  factor  of  5. 
Laser  velocimetry  measurements  by  Barata28  of  a  normal  impinging  jet  at  H/D  =  5  and 
V/V.  -  30  indicate  that  the  shear  stress  in  the  vortex  is  roughly  an  order  of  magnitude  less 
than  the  peak  values  in  the  free  jet.  In  contrast,  the  CVS  model  predicts  the  shear  stress  to 
be  small  in  the  ground  vortex  and  relatively  larger  in  the  free  jet.  Hence,  the  CVS  model 
displays  results  which  are  consistent  with  the  trends  seen  in  impinging  jet  experiments. 
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One  disturbing  characteristic  of  the  CVS  prediction  is  that  the  maximum  predicted 
shear  stress  occurs  just  downstream  of  the  jet  exit  on  the  upstream  side  of  the  jet. 
Experiments  indicate  the  shear  stress  is  small  near  the  nozzle  and  grows  in  the  downstream 
direction.  This  region  of  maximum  shear  stress  experiences  vortex  stretching,  because 
aW/dz  is  large  as  the  crossflow  is  deflected  around  the  normal  jet.  Furthermore,  the 
freestream  has  low  values  of  k  and  e.  This  combination  of  vortex  stretching  and  nearly 
laminar  flow  can  cause  the  model  to  display  ill-conditioned  behavior,  as  noted  previously. 
The  high  shear  stress  near  the  nozzle  exit  is  probably  due  to  this  problem  of  ill-conditioning. 
This  shortcoming  is,  in  some  respects,  an  inherent  disadvantage  of  all  multi-equation 
models  that  involve  ratios  of  turbulence  variables  which  become  singular  in  laminar  flow. 
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Figure  26.  Surface  centerline  pressure  distributions  computed  in  grid  refinement 
study. 

Grid  Dependence  -  The  grid  dependence  was  evaluated  for  this  case  by  comparing 
calculations  on  two  different  grids,  one  which  contained  50x35x28  points  in  the  x-  y-  and 
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z-directions  and  was  typical  of  the  grids  used  in  this  study,  and  one  which  contained 
94x46x44  grid  points.  In  the  finer  grid,  the  grid  spacing  was  approximately  half  of  that  in 
the  regular  grid,  in  the  free  jet,  the  impingement  region,  and  ground  vortex  region.  This 
should  reduce  the  truncation  by  a  factor  of  roughly  16  because  most  error  is  due  to  the 
inviscid  terms  which  are  discretized  with  fourth  order  accuracy.  The  calculations  were 
performed  with  the  VS  model,  although  the  specific  model  is  largely  irrelevant  in  a  grid 
refinement  study.  Figure  26  gives  the  centerline  surface  pressure  distribution  under  the 
impingement  point  and  ground  vortex.  The  biggest  difference  is  seen  in  the  peak  pressure, 
with  the  fine-grid  solution  giving  a  higher  pressure.  This  result  is  consistent  with  the 
reduced  "numerical  diffusion"  in  the  fine-grid  solution  of  the  free  jet.  The  position  of  the 
ground  vortex,  associated  with  the  low-pressure  region,  is  very  similar  for  the  calculations 
on  both  grids,  and,  consistent  with  the  reduced  losses  due  to  numerical  diffusion,  the 
ground  vortex  in  the  fine-grid  calculation  is  forward  of  that  on  the  regular  grid.  However, 
the  dependence  of  the  computed  vortex  location  upon  the  grid  is  significantly  less  than  the 
dependence  upon  the  turbulence  model,  and  turbulence  modeling  can  be  addressed  without 
significant  concern  regarding  truncation  errors. 

Case  II:  Normal  Jet,  H/D  =  3,  V/V_  =  10 

This  case  was  also  measured  by  Cimbala  et  al,29and  the  jet  configuration  is  similar  to 
Case  I,  except  for  the  velocity  ratio.  A  significant  part  of  the  experimental  study  focused  on 
vortex  unsteadiness,  and  the  core  vortex  height  fluctuates  ±30%  in  this  flow.  Figure  27  gives 
a  comparison  of  the  time-averaged  velocity  vectors  in  the  impinging  jet  and  ground  vortex 
region  computed  with  the  CVS  model  and  the  experimental  data.  The  location  of  the 
ground  vortex  core  is  measured  to  be  at  x/Dj  =  4.5  and  y/Dj  =  0.9.  The  computed  location  is 
at  x/Dj  =  5.25  and  y/D:  =  1.  The  computed  vortex  core  location  is  about  17%  too  far 
upstream,  compared  to  the  experimental  location. 


(a)  Measured  velocities 
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(b)  Predicted  velocities 


Figure  27 .  Comparison  of  center  plane  mean  velocity  vectors  for  Case  II  computed 

with  the  CVS  model  and  measured  by  Cimbala  et  al.29 


The  most  significant  other  difference  is  that  the  shear  layer  on  the  free  jet  is  thicker  in  the 
calculation  than  in  the  experiment.  There  are  between  5  and  9  grid  points  across  the  shear 
layer  at  the  plotted  locations,  which  should  be  sufficient  for  good  numerical  resolution  of 
the  mean  flow  with  the  present  differencing  scheme.  Hence,  the  excessive  spreading  is 
probably  caused  by  the  turbulence  model.  As  noted  previously,  the  standard  k-e  model 
overpredicts  the  spreading  rate  of  the  round  jet,  and  the  vortex  stretching  modification 
exacerbates  this  problem.  While  some  of  the  physics  of  this  issue  are  understood,  no 
suitable  model  correction  has  yet  been  proposed. 

Case  III:  Inclined  Impinging  Jet 

The  final  calculation  presented  is  of  the  primary  experimental  case,  which  is  an 
impinging  jet  inclined  45°  toward  the  crossflow.  The  jet  is  close  to  fully  developed 
turbulent  pipe  flow  with  a  jet  centerline  velocity  ratio  U«d/U.  =  5.27,  and  H/D  =  3.1.  The 
wind  tunnel  cross  section  is  76.2  cm  (30  in)  square  and  the  jet  exit  diameter  D-  =  2.87  cm 
The  nozzle  exit  centerline  speed  is  approximately  26  m/s.  The  developed  turbuient  jet  was 
used  in  this  case  because  an  important  application  of  inclined  jets  is  for  thrust  reversing;  the 
reversing  process  on  jet  engines  tends  to  generate  a  turbulent  jet. 
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(a)  Velocity  magnitude 


<b)  Velocity  vectors 


Figure  28.  Contours  of  velicity  magnitude,  velocity  vectors,  and  contours  of  the  shear 
stress  u  V  on  the  center  plane  for  Case  I II  predicted  with  the  CVS  model. 
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(c)  Shear  stress 


Figure  28.  Contours  of  velicity  magnitude,  velocity  vectors,  and  contours  of  the  shear 
b  stress  uV  on  the  center  plane  for  Case  III  predicted  with  the  CVS  model. 


Velocity  magnitude  contours,  velocity  vectors,  and  contours  of  the  shear  stress  u'v  ' 
predicted  with  the  CVS  model  are  given  in  Figure  28.  In  the  velocity  magnitude  contours 
and  velocity  vectors,  the  primary  features  of  the  flow  can  be  seen.  The  inclined  jet  has  no 
potential  core,  and  the  velocity  magnitude  decays  significantly  before  reaching  the  ground. 
The  majority  of  the  jet  flow  goes  forward  after  the  interaction  with  the  ground,  and  the 
forward  extent  of  the  ground  vortex,  defined  by  the  separation  line  between  the  crossflow 
and  the  ground  vortex,  is  about  x/Dj  -  -13.  The  contour  plot  of  u  v  must  be  interpreted 
correctly  because  the  impinging  jet  is  at  -45°  to  the  principal  grid  axis.  The  effectjye  shear 
stress"  which  transports  momentum  away  from  the  jet  does  not  appear  as  u'v',  but  is 
shared  by  u'u',  v'v',  and  u'v',  according  to  rules  of  tensor  rotation.  The  areas  of 
maximum  u'v',  are  the  windward  shear  layer  of  the  impinging  jet,  a  region  near  the 
downstream  edge  of  the  impingement  zone,  and  in  the  wall  jet  running  forward  from  the 
impingement  point.  As  in  Case  I,  there  is  negligible  shear  stress  in  the  vortex  core.  This 
plot  also  reveals  a  region  of  shear  stress  in  the  ground  vortex,  but  above  the  core. 

The  accuracy  of  the  model  predictions  is  best  assessed  in  a  plot  of  the  centerline  surface 
pressure  coefficient,  given  in  Fig.  29.  For  this  flow,  the  STD  model  underpredicts  the 
forward  extent  of  the  ground  vortex,  defined  by  the  location  of  minimum  pressure,  by 
roughly  30%,  while  the  CVS  model  predicts  the  position  of  the  ground  vortex  core  to  be 
roughly  1 8%  upstream  of  the  measured  value.  The  VS  model  causes  the  ground  vortex  to 
move  backward  slightly  from  the  location  predicted  with  the  STD  model,  further  degrading 
the  accuracy.  This  result  clearly  indicates  that  the  apparent  accuracy  of  the  VS  model  in 
Case  I  is  a  result  of  luck,  rather  than  fundamentally  correct  turbulence  modeling. 


(  X  ~  X  nozzle  ) 


Figure  29.  Surface  centerline  static  pressures  for  Case  III. 


There  are  significant  differences  in  the  predicted  maximum  wall  pressures  under  the 
impingement  region,  and  the  specific  causes  of  the  differences  are  not  known.  Several 
factors  affect  this  pressure.  The  vortex  stretching  model  may  overpredict  the  turbulent 
stresses  in  the  round  jet  above  the  impingement  point,  which  affects  the  jet  deceleration  and 
deflection  by  the  crossflow.  The  maximum  static  wall  pressure  is  not  the  maximum  jet 
stagnation  pressure,  because  the  stagnation  streamline  does  not  correspond  to  the  jet 
centerline  in  an  inclined  jet.  Hence,  differences  in  the  predicted  turbulence  may  generate 
significant  differences  in  the  maximum  surface  pressure,  and  the  relationship  between  the 
pressure  and  predicted  turbulence  is  not  obvious. 

Interpretation  of  Reynolds-Averaged  Results 

The  above  results  provide  a  self-consistent  description  of  the  deficiencies  of  the 
standard  k-£  model,  and  of  the  improvements  in  accuracy  given  by  modifications  to  the 
standard  model  to  account  for  effects  of  curvature  and  vortex  stretching. 
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The  STD  model  predicted  the  ground  vortex  to  be  50%  forward  of  its  measured 
location  in  Case  I  and  30%  downstream  in  Case  III.  This  80%  discrepancy  in  the  predicted 
relative  ground  vortex  locations  indicates  that  powerful  turbulence-modifying  phenomena 
occur  in  the  impingement  region,  which  depend  significantly  upon  the  angle  of  the  jet.  The 
VS  modification  produced  changes  in  Cases  I  and  III  that  were  of  the  same  sign;  the  same 
was  true  for  the  C  modification.  However,  the  CVS  modification  causes  the  ground  vortex 
to  move  downstream  in  Case  I  and  upstream  in  Case  III,  which  indicates  that  the  predicted 
relative  effects  on  turbulence  of  vortex  stretching  and  curvature  are  different  in  the  two 
cases. 

The  wall  jet  upstream  of  the  impingement  point  spreads  out  more  rapidly  in  the 
normal  impinging  jet  than  in  the  inclined  impinging  jet.  Hence,  3W/3z,  appropriately 
normalized,  is  larger  for  the  normal  impinging  jet,  and  vortex  stretching  effects  in  the  wall 
jet  are  stronger  in  Case  I  than  in  Case  III.  Differences  in  streamline  curvature  in  Cases  I  and 
III  are  appreciable  in  the  impingement  region,  and  should  produce  differences  in  the 
modeled  (and  physical)  turbulence  there.  However,  the  region  of  strong  curvature  in  the 
impingement  zone  is  small  relative  to  the  extent  of  the  radial-like  wall  jet.  Furthermore, 
turbulence  in  the  vortex  core  is  suppressed,  and  the  predicted  turbulent  shear  stress  in  the 
vortex  cores  approaches  zero  in  Cases  I  and  III.  There  is  no  comparable  limit  to  turbulence 
amplification  in  the  wall  jet.  Hence,  the  significant  differences  seen  in  the  accuracy  of  the 
STD  model  for  the  different  impinging  jet  flows  is  apparently  due  to  that  model's  inability 
to  predict  the  effects  of  vortex  stretching  upstream  of  the  impingement  point. 

The  CVS  model  gave  a  significant  improvement  over  the  STD  model  in  Cases  I  and  II, 
and  consistently  predicted  the  ground  vortex  location  in  all  three  cases  considered  here  to 
be  13%  - 18%  upstream  of  the  measured  location.  Hence,  the  errors  in  solutions  predicted 
by  the  CVS  model  are  one-half  to  one-third  of  the  errors  given  by  the  STD  model,  and  they 
are  consistent  across  the  cases  examined. 

The  errors  given  by  the  CVS  model  are  not,  however,  small  enough  to  be  disregarded, 
and  some  means  of  reducing  these  errors  is  desired.  In  the  use  of  most  turbulence  models, 
the  model  coefficients  can  be  adjusted  to  tune  the  model  predictions  for  a  certain  class  of 
flows.  The  calculations  presented  above  were  obtained  with  model  coefficients  tuned  for 
benchmark  flows  which  were  totally  distinct  from  impinging  jet  flows.  It  is  legitimate  to 
change  these  coefficients,  provided  the  coefficients  can  be  held  constant  for  the  entire  class 
of  flows  of  interest.  An  attempt  was  made  to  optimize  the  coefficients  Cvs  and  Cvs2  for  flows 
including  the  wall  jet  and  impinging  jet  flows.  The  coefficients  were  changed  to  increase 
effects  of  near-wall  vortex  stretching  and  reduce  its  effects  in  regions  far  from  a  wall.  One 
objective  of  this  test  was  to  reduce  the  very  high  stresses  seen  near  the  nozzle  exit  in  Case  I, 
However,  the  results  of  this  test  were  negative;  it  was  not  possible  to  improve  the  predicted 
ground  vortex  location  with  a  pair  of  coefficients  that  would  also  predict  the  radial  wall  jet 
accurately.  The  inability  to  obtain  improved  accuracy  in  a  range  of  flows  by  adjusting 
model  coefficients  suggests  that  the  underlying  physical  mechanisms  must  be  represented 
more  accurately  in  the  turbulence  model. 
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7.  CONCLUSION 


The  present  study  addresses  the  accuracy  of  the  k— e  turbulence  model  for  inclined 
impinging  jet  flows  typical  of  thrust  reverser  application  in  close  ground  proximity, 
particularly  those  which  generate  a  complex  recirculating  flow  field. 

The  thesis  of  the  present  report  is  that  streamline  curvature  and  vortex  stretching  are 
the  two  most  significant  influences  on  turbulence  in  these  flows  which  are  not  treated 
adequately  in  the  "standard"  k—  t  model.  Streamline  curvature  has  been  recognized  as  a 
critical  turbulence  modifier  for  many  years.  The  significance  of  vortex  stretching  in 
impinging  jet  flows  was  recognized  in  the  large  eddy  simulations  performed  in  this  study. 
As  vortical  structures  flow  outward  from  the  impingement  point  in  a  radial-like  wall  jet,  the 
structures  contract  and  intensify  due  to  conservation  laws  for  mass  and  angular  momentum. 
Instantaneous  and  time-averaged  results  from  the  large  eddy  simulations  reveal  that  the 
shear  stress  is  amplified  by  the  vortex  stretching  process. 

Modifications  which  account  for  curvature  and  vortex  stretching  were  added  to  the 
standard  k—  e  model,  and  the  new  coefficients  in  the  resulting  modified  model,  denoted  the 
CVS  model,  were  determined  via  numerical  optimization  for  benchmark  jet  flows  unrelated 
to  the  impinging  jet  flows  of  primary  concern. 

In  the  three  cases  considered,  two  impinging  jets  issuing  from  a  nozzle  normal  to  the 
ground  plane  and  one  issuing  from  a  nozzle  inclined  45*  into  the  crossflow,  the  CVS  model 
consistently  improved  the  accuracy  for  the  computed  location  of  the  ground  vortex.  Errors 
in  the  computed  vortex  location,  relative  to  the  distance  between  the  impingement  point 
and  the  ground  vortex  core,  were  as  large  as  50%  in  calculations  with  the  standard 
k-e  model  and  could  be  of  either  sign  (overprediction  or  underprediction  of  the  ground 
vortex's  forward  extent).  The  CVS  model  reduced  these  errors  to  the  range  of  13%-18%,  and 
it  consistently  overpredicted  the  upstream  penetration  of  the  ground  vortex. 

Recommended  Future  Work 

There  are  a  few  significant  shortcomings  in  the  turbulence  model  at  this  time.  The  first 
is  a  flaw  in  the  STD  model:  it  cannot  accurately  predict  the  "simplest"  elements  of  thrust 
reverser  flow  fields.  It  does  not  properly  predict  a  planar  wall  jet  on  a  flat  surface  or  a 
round  jet  if  the  model  coefficients  have  been  adjusted,  as  is  normally  done,  to  predict  the 
widest  range  of  generic  turbulent  flows.  The  standard  k-f  model  simply  does  not  have  the 
capacity  to  respond  to  the  physically  relevant  parameters  such  as  wall  proximity,  in  the  wall 
jet,  or  the  change  of  dominant  vortical  structure,  as  in  the  round  jet.  These  errors  are  of  the 
order  of  10%  to  20%,  which  is  large  enough  to  make  a  significant  contribution  to  the  total 
error,  or  mask  other  errors,  in  the  calculation  of  a  complex  impinging  jet  flow. 

An  additional  problem  is  that  the  k-£  model,  as  with  many  multi-equation  models,  is 
singular  in  laminar  flow.  Although  the  ratio  k/ e  may  be  well-defined  physically,  even 
small  truncation  errors  can  generate  non-physical  values  of  k  and  t  which  will  cause  the 
modeled  solution  to  behave  erratically.  If  the  objective  of  developing  good  turbulence 
models  for  complex  flows  is  to  be  achieved,  these  models  must  be  well-posed  for  integration 
by  numerical  schemes  of  finite  accuracy. 
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Modeling  of  the  vortex  stretching  process  may  also  benefit  from  further  work.  The 
present  study  has  identified  vortex  stretching  as  a  critical  physical  aspect  of  impinging  jet 
flows  and  has  employed  a  relatively  simple  model  for  this  process.  However,  it  has  not 
been  possible  to  obtain  a  good  understanding  of  the  physics  of  vortex  stretching,  nor  has  it 
been  possible  to  develop  a  very  good  model  of  this  physics. 

It  is  not  difficult  to  make  general  recommendations  for  research  work  that  will  address 
the  above  problems:  employ  turbulence  simulations  and  experiments  symbiotically  to 
identify  the  dominant  turbulent  physics,  discriminate  between  that  which  can  be  modeled 
simply  and  that  which  is  complex,  and  construct  a  mathematically  well-posed  model  that 
resolves  that  dominant,  complex  physics.  Unfortunately,  decades  of  research  have  not  yet 
achieved  these  goals. 

Specific  recommendations  to  improve  k-«  modeling  for  thrust  reverser  flows  are 
somewhat  more  difficult  to  put  forth.  Clearly,  the  effects  of  vortex  stretching  on  turbulence 
should  be  the  topic  of  a  basic  research  effort.  The  roles  of  the  inherent  vortex  structure,  the 
modification  of  that  structure  by  stretching,  and  the  subsequent  effect  on  turbulent  mixing 
must  be  known.  These  results  should  give  a  strong  indication  of  the  correct  direction  for 
future  modeling  work.  It  is  believed  that  the  k-e  model  carries  too  little  information  about 
turbulence  to  permit  it  to  make  good  predictions  of  the  wide  range  of  flow  phenomena 
which  occur  in  thrust  reverser  flow  fields.  This  does  not  imply  that  the  next  conventional 
step  in  the  turbulence  modeling  hierarchy,  the  Reynolds-stress  transport  model  (RSTM),  is 
the  appropriate  model  to  pursue.  The  RSTM's  accuracy  for  the  plane  wall  jet,  for  example, 
is  worse  than  that  of  the  standard  k—  t  model.  Instead,  it  is  believed  that  a  model  which 
contains  information  about  the  dominant  vortical  structure  will  be  required  to  achieve  a 
significant  improvement  over  existing  models. 
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